1. 概要
単位的可換環 $R$ 上の数列 $A = (A_0,A_1,A_2,\ldots)$ が,定数 $c_1, c_2, \ldots, c_d$ について
$$ A _ i = c_1 A _ {i-1} + c _ 2A _ {i-2} + \cdots + c _ dA _ {i-d}\qquad(i\geq d) $$
を満たすとします.このような数列を線形漸化的数列といい,次の講座でその性質を詳しく整理しました.
本記事では上の講座の続きとして,線形漸化式と初期値 $A_0, A_1, \ldots, A_{d-1}$ が与えられたときに,数列の第 $K$ 項 $A_K$ を求める問題をより詳しく扱います.
この問題を解く方法のひとつは行列累乗を用いるもので,単純な方法では計算量が $\mathrm{O}(d ^ 3 \log K)$ 時間となります.本記事では,多項式を用いてより高速に第 $K$ 項を求める方法を説明します.より具体的には,次数 $d$ 程度の多項式積の計算量を $M(d)$ としたとき,$\mathrm{O}(M(d)\log K)$ 時間で $A_K$ を求める方法として,Fiduccia のアルゴリズムと Bostan–Mori のアルゴリズムを扱います.
2. 前提知識
以下の講座を理解していることを前提とします.
また,多項式の乗算・剰余算を扱います.高速な実装では FFT を用いた多項式乗算が必要になりますが,アルゴリズムの考え方を理解するだけであれば,素朴な多項式演算を考えれば十分です.
3. 問題設定
本記事で扱うのは次の問題です.
数列 $A=(A_0,A_1,A_2,\ldots)$ が,定数係数線形漸化式 $$ A _ i = c _ 1A _ {i-1} + c _ 2A _ {i-2} + \cdots + c _ dA _ {i-d}\qquad(i\geq d) $$ を満たしているとします.漸化式の係数 $c_1, c_2, \ldots, c_d$,初期値 $A_0, A_1, \ldots, A_{d-1}$ および非負整数 $K$ が与えられたとき,数列 $A$ の第 $K$ 項 $A_K$ を求めてください.
ただし,次の約束をします.
- 単位的可換環 $R$ をひとつ固定し,数列の要素や漸化式の係数はすべて $R$ の元とします.
- $R$ の環演算(加算・減算・乗算)は $\mathrm{O}(1)$ 時間で行えるものとします.
代数の用語に不慣れな場合,$R=\mathbb{Z}/n\mathbb{Z}$ または $\mathbb{F}_p$ と考えてください.
また,計算量解析において,次数 $d$ の多項式積を計算する時間を $M(d)$ と書きます.例えば素朴な多項式乗算では
$$ M(d) = \mathrm{O}(d ^ 2) $$
です.競技プログラミングにおいて頻出する $R$($\mathbb{Z}/n\mathbb{Z}$,$\mathbb{F}_p$,$\mathbb{C}$ など)では,適切な計算量モデルのもとで,$M(d)=\mathrm{O}(d\log d)$ と考えることができます.
以降では,第 $K$ 項 $A_K$ を $\mathrm{O}(M(d)\log K)$ 時間で求める方法を説明します.
4. Fiduccia のアルゴリズム
本節では,線形漸化的数列の第 $K$ 項を求める古典的な方法である Fiduccia のアルゴリズムを説明します.
4.1. 線形漸化式と多項式剰余
数列 $A$ が,定数係数線形漸化式
$$ A _ i = c _ 1A _ {i-1} + c _ 2A _ {i-2} + \cdots + c _ dA _ {i-d}\qquad(i\geq d) $$
を満たすとします.この式をそのまま「左辺を右辺に置き換える操作が行える」と解釈すると,$A_i$ をより添字が小さい項の線形結合として表すことができます.このことを繰り返すと,$A_i$ を $A_0, A_1, \ldots, A_{d-1}$ の線形結合として表すことができます.
具体例
例えば $A$ が
$$ A_i = 2 A _ {i-1} + 3 A _ {i-2} \qquad (i\geq 2) $$
を満たすとします.$A _ 4$ より添字が小さい項の線形結合で表してみましょう.まず漸化式で $i=4$ とした式から
$$ A _ 4 = 2 A _ 3 + 3 A _ 2 $$
となります.次に漸化式で $i=3$ とした式 $A _ 3 = 2 A _ 2 + 3 A _ 1$ を用いると,
$$ A _ 4 = 2 (2 A _ 2 + 3 A _ 1) + 3 A _ 2 = 7 A _ 2 + 6 A _ 1 $$
となります.さらに漸化式で $i=2$ とした式 $A _ 2 = 2 A _ 1 + 3 A _ 0$ を用いると,
$$ A _ 4 = 7 (2 A _ 1 + 3 A _ 0) + 6 A _ 1 = 20 A _ 1 + 21 A _ 0 $$
となります.これで $A_4$ が $A _ 0, A _ 1$ の線形結合として表すことができました.
このように $A_i$ を,$A _ 0, A _ 1, \ldots, A _ {d-1}$ の線形結合で表すことができれば,第 $i$ 項を初期値から求めることができます.Fiduccia のアルゴリズムでは,この線形結合の係数を多項式剰余として求めます.
多項式剰余との比較
多項式 $f(x), g(x), \Gamma(x)$ について,$f(x)-g(x)=\Gamma(x)h(x)$ であるような多項式 $h(x)$ が存在するときに $f(x)\equiv g(x)\pmod{\Gamma(x)}$ と表すことにします.上の具体例の計算過程は,次のような多項式剰余を求める計算過程と類似していることが分かります.
$$\Gamma(x) = x ^ 2 - (2 x + 3)$$
とします.$x ^ 4$ を $\Gamma(x)$ で割った余りは次のように求められます.
$$x ^ 2 \equiv 2 x + 3\pmod{\Gamma(x)} \tag{$*$}$$
であることに注意します.
まず $(*)$ の $x ^ 2$ 倍から
$$ x ^ 4 \equiv 2 x ^ 3 + 3 x ^ 2 \pmod{\Gamma(x)} $$
となります.次に $(*)$ の $x$ 倍から分かる式 $x ^ 3 \equiv 2 x ^ 2 + 3 x\pmod{\Gamma(x)}$ を用いると,
$$ x ^ 4 \equiv 2 (2 x ^ 2 + 3 x) + 3 x ^ 2 = 7 x ^ 2 + 6 x \pmod{\Gamma(x)} $$
となります.さらに $(*)$ を用いると,
$$ x ^ 4 \equiv 7 (2 x + 3) + 6 x \equiv 20 x + 21 \pmod{\Gamma(x)} $$
となります.これで $x ^ 4 \equiv r _ 0 + r _ 1 x \pmod{\Gamma(x)}$ を満たす $r _ 0, r _ 1$ (言い換えれば,$x ^ 4$ の $\Gamma(x)$ による剰余)を求めることができました.
定理
上の具体例で見た通り,
- 線形漸化的数列の $A _ i$ を $A _ 0, A _ 1, \ldots, A _ {d-1}$ の線形結合で表す問題
- $x ^ i$ の $\Gamma(x)$ による剰余 $r _ 0 + r _ 1 x + \cdots + r _ {d-1} x ^ {d-1}$ を求める問題
はほとんど同じ計算により解くことができます.ここから以下の定理が得られます.
定数係数線形漸化式
$$ A _ i = c _ 1 A _ {i-1} + c _ 2 A _ {i-2} + \cdots + c _ d A _ {i-d} \qquad (i\geq d) $$に対して,多項式 $$ \Gamma(x)=x ^ d - c _ 1 x ^ {d-1} - c _ 2 x ^ {d-2} - \cdots - c _ d $$ を特性多項式(characteristic polynomial)という.
数列 $A=(A_0,A_1,A_2,\ldots)$ を位数 $d$ の定数係数線形漸化式
$$ A _ i = c _ 1 A _ {i-1} + c _ 2 A _ {i-2} + \cdots + c _ d A _ {i-d} \qquad (i\geq d) $$を満たす線形漸化的数列とし,
$$ \Gamma(x)=x ^ d - c _ 1 x ^ {d-1} - c _ 2 x ^ {d-2} - \cdots - c _ d $$をその特性多項式とする.非負整数 $i$ に対して $x ^ i$ の $\Gamma(x)$ による剰余を
$$ r _ 0 + r _ 1 x + \cdots + r _ {d-1} x ^ {d-1} $$とするとき,
$$ A _ i = r _ 0 A _ 0 + r _ 1 A _ 1 + \cdots + r _ {d-1} A _ {d-1} $$が成り立つ.
上の具体例で見たように,この定理は,$x ^ d \equiv c_1x ^ {d-1} + \cdots + c_d \pmod{\Gamma(x)}$ という置き換えと,$A _ i = c_1A _ {i-1} + \cdots + c _ d A _ {i-d}$ という置き換えが同じ係数を持つことから従います.$\blacksquare$
数列 $A$ を漸化式 $$ A _ {i} = A _ {i-1} + A_{i-2} + A_{i-3} \qquad (i\geq 3) $$ および初期値 $A _ 0 = 3, A _ 1 = 2, A _ 2 = 1$ で定まる数列とします.
- 漸化式を用いて $A_3, A_4, A_5$ を求めてください.
- $A_5$ を $A_0, A_1, A_2$ の線形結合として表し,そこから $A_5$ を求めてください.
解答
- $A _ 3 = 6, A _ 4 = 9, A _ 5 = 16$.
- $A _ 5 = 4 A _ 2 + 3 A _ 1 + 2 A _ 0$.ここから $A _ 5 = 4\cdot 1 + 3\cdot 2 + 2\cdot 3 = 16$.
4.2. Fiduccia のアルゴリズム
前節の定理より,第 $K$ 項 $A_K$ を求めるには,$x ^ K$ の特性多項式 $\Gamma(x)$ による剰余を求めればよいことが分かりました.
$x ^ K \bmod \Gamma(x)$ は,繰り返し二乗法により計算できます.これを用いると,次の擬似コードによって $A _ K$ を求めることができます.
function Fiduccia(A, Gamma(x), K):
# A[0], A[1], ..., A[d-1]: initial values
# Gamma(x): characteristic polynomial of degree d
# returns A_K
function ModPow(K):
# returns x^K mod Gamma(x)
if K < d:
return x^K
R(x) := ModPow(floor(K / 2))
R(x) := R(x)^2 mod Gamma(x)
if K is odd:
R(x) := x R(x) mod Gamma(x)
return R(x)
R(x) := ModPow(K)
ans := 0
for i = 0, 1, ..., d - 1:
ans := ans + ([x^i]R(x)) * A[i]
return ans
ここで,次数 $d$ 未満の多項式どうしの積を $\Gamma(x)$ で割った余りを求める操作は,素朴な方法で $\mathrm{O}(d ^ 2)$ 時間で行えます.また本記事では解説を省略しますが,畳み込みと形式的べき級数逆元を用いて $\mathrm{O}(M(d))$ 時間で行えることが知られています.
したがって,高速な多項式剰余算を用いることで,Fiduccia のアルゴリズムにより,線形漸化的数列の第 $K$ 項 $A_K$ を $\mathrm{O}(M(d)\log K)$ 時間で求めることができます.
$2$ 種の繰り返し二乗法による違い
繰り返し二乗法のアルゴリズムとしては,次のアルゴリズムも有名です.
function ModPow(K):
# returns x^K mod Gamma(x)
R(x) := 1
B(x) := x
while K > 0:
if K is odd:
R(x) := R(x) B(x) mod Gamma(x)
B(x) := B(x)^2 mod Gamma(x)
K := floor(K / 2)
return R(x)
いずれの実装でも,$\mathrm{O}(M(d) \log K)$ 時間で $x ^ K \bmod \Gamma(x)$ を求めることができます.ただし,計算量には定数倍の違いが生じます.
本注釈の擬似コードのアルゴリズムでは,$K$ の各ビットごとに,$\Gamma(x)$ を法とする多項式乗算が最大 $2$ 回生じます.本文中の擬似コードでも,見かけ上は各段で最大 $2$ 回の多項式乗算が生じます.
しかし本文中の擬似コードはより詳細に見れば,$2$ 回のうち $1$ 回は
R(x) := x R(x) mod Gamma(x)
です.同じ「$\Gamma(x)$ を法とする多項式乗算」でも,この場合には乗算する片方の多項式が $x$ に限られています.この計算は係数をひとつずらしたあと,$x ^ d$ の項を
$$ x ^ d \equiv c _ 1 x ^ {d-1} + c _ 2 x ^ {d-2} + \cdots + c _ d \pmod{\Gamma(x)} $$
で置き換えればよく,$\mathrm{O}(d)$ 時間で計算できます.そのため Fiduccia のアルゴリズムでは,本注釈の擬似コードよりも,本文中の擬似コードの形を用いる方が定数倍の面で優れています.
呼称について
4.3. 実装例
Typical DP Contest「T - フィボナッチ」 への提出です.
この実装では,多項式乗算・剰余算を素朴な方法により $\mathrm{O}(d ^ 2)$ 時間で計算しています.漸化式の位数を $d$,求める添字を $K$ とすると,解法全体の計算量は $\mathrm{O}(d ^ 2\log K)$ 時間です(問題文の記号では,$\mathrm{O}(K ^ 2\log N)$ 時間です).
5. Bostan–Mori のアルゴリズム
本節では,線形漸化的数列の第 $K$ 項を求めるもうひとつの方法として,Bostan–Mori のアルゴリズムを解説します.
5.1. 線形漸化式と形式的べき級数(復習)
線形漸化的数列 で解説したように,線形漸化的数列 $A = (A_0, A_1, A_2, \ldots)$ の母関数
$$ A(x) = A _ 0 + A _ 1 x + A _ 2 x ^ 2 + \cdots $$
は有理式で表すことができます.数列 $A$ が定数係数線形漸化式
$$ A _ i = c _ 1A _ {i-1} + c _ 2A _ {i-2} + \cdots + c _ dA _ {i-d}\qquad(i\geq d) $$
を満たすとします.このとき,母関数 $A(x)$ は $A(x)=\dfrac{P(x)}{Q(x)}$ と表せます.ここで
$$ \begin{aligned} Q(x) &= 1 - c _ 1 x - c _ 2 x ^ 2 - \cdots - c _ d x ^ d, \\ P(x) &= A(x)Q(x) \bmod x ^ d \end{aligned} $$
です.第 $K$ 項 $A_K$ は $A_K=[x ^ K]\dfrac{P(x)}{Q(x)}$ と表せます.Bostan–Mori のアルゴリズムは,この有理式の第 $K$ 係数を高速に求めるアルゴリズムです.
なお $Q(x)$ は定数係数線形漸化式の係数を並べてできる多項式ですが,特性多項式とは係数の並び順が逆順であることに注意してください.
5.2. 形式的べき級数の偶数次部分,奇数次部分
形式的べき級数
$$ F(x) = \sum _ {i = 0} ^ {\infty} F _ i x ^ i = F _ 0 + F _ 1 x + F _ 2 x ^ 2 + \cdots $$
に対して,その偶数次部分,奇数次部分を
$$ \begin{aligned} F _ {\mathrm{e}}(x) &= \sum _ {i = 0} ^ {\infty} F _ {2i} x ^ i = F _ 0 + F _ 2 x + F _ 4 x ^ 2 + \cdots, \\ F _ {\mathrm{o}}(x) &= \sum _ {i = 0} ^ {\infty} F _ {2i+1} x ^ i = F _ 1 + F _ 3 x + F _ 5 x ^ 2 + \cdots. \end{aligned} $$
により定めます.$F(x) = F _ {\mathrm{e}}(x ^ 2) + x F _ {\mathrm{o}}(x ^ 2)$ が成り立ちます.
Bostan–Mori のアルゴリズムでは,次の補題を用います.
$Q(x)$ を形式的べき級数とするとき,$Q(x)Q(-x)$ の奇数次部分は $0$ である.つまりある形式的べき級数 $V(x)$ が存在して,
$$ Q(x)Q(-x)=V(x ^ 2) $$が成り立つ.
非負整数 $i$ に対して,$[x ^ i]Q(x)$ を $Q_i$ と書くことにします.非負整数 $k$ に対して
$$ [x ^ k] Q(x) Q(-x) = \sum _ {i + j = k} \left([x ^ i]Q(x)\right)\cdot \left([x ^ j]Q(-x)\right) = \sum _ {i + j = k} Q _ i\cdot (-1) ^ j Q_j $$
が成り立ちます.さらに $k$ が奇数であるとき,$i+j=k$ を満たす非負整数の組 $(i,j)$ は次のようなものです.
- $a + b = k$ を満たす偶数 $a$,奇数 $b$ に対して $(i,j) = (a,b)$.
- $a + b = k$ を満たす偶数 $a$,奇数 $b$ に対して $(i,j) = (b,a)$.
同じ $(a,b)$ に対するこれら $2$ 通りの $(i,j)$ の寄与をまとめると,
$$ Q _ a\cdot (-1) ^ b Q_b + Q _ b \cdot (-1) ^ a Q _ a = (-Q _ a Q _ b) + Q _ bQ_a = 0 $$
となるので,$[x ^ k] Q(x)Q(-x)=0$ が成り立ちます.$\blacksquare$
証明方法に関する補足
5.3. Bostan–Mori のアルゴリズム
定数係数線形漸化式を満たす数列の第 $K$ 項を求める問題は,
- $d-1$ 次以下の多項式 $P(x)$
- $d$ 次以下の多項式 $Q(x)$(ただし $[x ^ 0]Q(x)=1$)
に対して $[x ^ K]\dfrac{P(x)}{Q(x)}$ を求める問題に帰着できました.Bostan–Mori のアルゴリズムはこの問題を,次のように再帰的に解くものです.
まず,
$$ \frac{P(x)}{Q(x)}= \frac{P(x)Q(-x)}{Q(x)Q(-x)} $$
と変形します.【補題 5】 より,分母はある $d$ 次以下の多項式 $V(x)$ により $V(x ^ 2)$ と書けます.分子を $U(x) = P(x)Q(-x)$ とおきます.すると
$$ \frac{P(x)}{Q(x)} = \frac{U _ {\mathrm{e}}(x ^ 2)}{V(x ^ 2)} + x\frac{U _ {\mathrm{o}}(x ^ 2)}{V(x ^ 2)} $$
が成り立ちます.この式は,$\dfrac{P(x)}{Q(x)}$ の偶数次部分,奇数次部分への分解に対応します.つまり,$\dfrac{P(x)}{Q(x)}$ の偶数次部分は $\dfrac{U _ {\mathrm{e}}(x)}{V(x)}$,奇数次部分は $\dfrac{U _ {\mathrm{o}}(x)}{V(x)}$ となります.
この式から,$K$ の偶奇に応じて次が分かります.
$$ \begin{aligned} [x ^ K]\frac{P(x)}{Q(x)} &= [x ^ {K/2}]\frac{U _ {\mathrm{e}}(x)}{V(x)} & (K\equiv 0\pmod{2}), \\ [x ^ K]\frac{P(x)}{Q(x)} &= [x ^ {(K-1)/2}]\frac{U _ {\mathrm{o}}(x)}{V(x)} & (K\equiv 1\pmod{2}). \end{aligned} $$
ここで右辺の分子に現れる多項式 $U _ {\mathrm{e}}(x), U _ {\mathrm{o}}(x)$ は $d-1$ 次以下の多項式です.また右辺の分母に現れる多項式 $V(x)$ は $d$ 次以下の多項式で,$[x ^ 0]V(x)=1$ を満たします.これで問題の形を保ったまま,求める係数の添字 $K$ を半分以下にすることができました.
この操作を繰り返すと,最終的には $K=0$ の場合,つまり $[x ^ 0]\dfrac{P(x)}{Q(x)}$ を求める問題に帰着されますが,$[x ^ 0]Q(x) = 1$ よりこの問題の答えは $[x ^ 0]P(x)$ です.
まとめると,以下の擬似コードにより $A_K$ を求めることができます.
function Bostan_Mori(P(x), Q(x), K):
# returns A_K = [x ^ K]P(x)/Q(x)
# assume [x^0]Q(x)=1
while K > 0:
U(x) := P(x) Q(-x)
W(x) := Q(x) Q(-x)
if K is even:
P(x) := U_e(x)
else:
P(x) := U_o(x)
Q(x) := W_e(x)
K := floor(K / 2)
return [x^0]P(x)
これが Bostan–Mori のアルゴリズムです.
イテレーションごとに多項式乗算 $P(x)Q(-x)$,$Q(x)Q(-x)$ を計算するため,計算量は $\mathrm{O}(M(d)\log K)$ 時間となります.
数列 $F = (F _ 0, F _ 1, F _ 2, \ldots) = (0, 1, 1, 2, 3, 5, 8, 13, 21, \ldots)$ を Fibonacci 数列とすると,$F$ の母関数は $F(x) = \dfrac{x}{1 - x - x ^ 2}$ となるのでした.
- $F _ {\mathrm{o}}(x) = \dfrac{U(x)}{V(x)}$ を満たす多項式 $U(x), V(x)$ を求めてください.
- 上で求めた $U(x), V(x)$ について $\dfrac{U(x)}{V(x)}$ について低次の係数を計算することで,$\dfrac{U(x)}{V(x)} = F _ 1 + F _ 3 x + F _ 5 x ^ 2 + F _ 7 x ^ 3 + \cdots$ となっていることを確かめてください.
解答
- $\dfrac{x}{1 - x - x ^ 2} = \dfrac{x (1 + x - x ^ 2)}{(1 - x - x ^ 2)(1 + x - x ^ 2)} = \dfrac{x + x ^ 2 - x ^ 3}{1 - 3 x ^ 2 + x ^ 4}$ より,奇数次部分は $\dfrac{1 - x}{1 - 3x + x ^ 2}$ となります.つまり $U(x) = 1 - x$,$V(x) = 1 - 3x + x ^ 2$ が条件を満たします.
- $A(x) = \dfrac{1 - x}{1 - 3x + x ^ 2}$ の係数 $A _ i$ が漸化式 $A _ i = 3 A _ {i-1} - A _ {i-2}$ を満たすことから,$\dfrac{1 - x}{1 - 3x + x ^ 2} = 1 + 2x + 5 x ^ 2 + 13 x ^ 3 + \cdots$ であることが確かめられます.
5.4. 実装例
Typical DP Contest「T - フィボナッチ」 への提出です.
この実装では,Bostan–Mori のアルゴリズムを用いています.多項式乗算は素朴な方法により計算しているため,漸化式の位数を $d$,求める添字を $K$ とすると,解法全体の計算量は $\mathrm{O}(d ^ 2\log K)$ 時間です(問題文の記号では,$\mathrm{O}(K ^ 2\log N)$ 時間です).
もうひとつの実装例を挙げます.これは Library Checker「Kth term of Linearly Recurrent Sequence」 への提出です.多項式乗算を FFT による畳み込みで行うことで,計算量は $\mathrm{O}(d\log d\log K)$ 時間になります.畳み込みの実装には AtCoder Library を用いています.
6. 関連問題
関連問題として挙げた問題は,一部 線形漸化的数列 と重複しています.
- https://atcoder.jp/contests/tdpc/tasks/tdpc_fibonacci
- https://yukicoder.me/problems/no/891
- https://atcoder.jp/contests/math-and-algorithm/tasks/math_and_algorithm_av
- https://atcoder.jp/contests/abc009/tasks/abc009_4
- https://judge.yosupo.jp/problem/kth_term_of_linearly_recurrent_sequence
- https://yukicoder.me/problems/no/215
- https://atcoder.jp/contests/abc178/tasks/abc178_d
- https://atcoder.jp/contests/s8pc-3/tasks/s8pc_3_g
- https://atcoder.jp/contests/ajo2025-final/tasks/ajo2025_final_b
- https://atcoder.jp/contests/abc436/tasks/abc436_g
- https://atcoder.jp/contests/abc300/tasks/abc300_h
- https://atcoder.jp/contests/fps-24/tasks/fps_24_t
問題について
1, 2, 3, 4, 5 は漸化式が与えられていて第 $K$ を求める問題です.このうち 1, 2, 3, 4 は制約が小さく,$O(d ^ 2)$ 時間での多項式演算で解けます.5 は FFT を用いた高速化が想定されています.
6, 7, 8, 9, 10 は比較的簡単に漸化式(あるいは母関数の有理式表示)が得られる問題です.なお,問題によっては負の二項定理などを用いてより高速に解くことも可能です.
11 は Bostan–Mori のアルゴリズムの応用です.12 は母関数の有理式表示を得るのが少し難しい問題です.
7. まとめ
本記事では,定数係数線形漸化式と初期値が与えられたときに,数列の第 $K$ 項を高速に求める方法として,Fiduccia のアルゴリズムと Bostan–Mori のアルゴリズムを解説しました.
Fiduccia のアルゴリズムは,線形漸化的数列の第 $K$ 項計算を,特性多項式 $\Gamma(x)$ による剰余 $x ^ K \bmod \Gamma(x)$ の計算と結びつける方法です.この考え方は,線形漸化的数列と,同伴行列・行列累乗・行列の特性多項式などの関係を理解する上でも重要です.
一方,Bostan–Mori のアルゴリズムは,線形漸化的数列の第 $K$ 項計算を母関数の係数 $[x ^ K]\dfrac{P(x)}{Q(x)}$ と考えて計算する方法です.
どちらの方法も,多項式乗算の計算量を $M(d)$ とすると,$\mathrm{O}(M(d)\log K)$ 時間で第 $K$ 項を求めることができます.特に,競技プログラミングで頻出する $R = \mathbb{Z} / n\mathbb{Z}$ における計算では $M(d)=\mathrm{O}(d\log d)$ と考えることができ,計算量は $\mathrm{O}(d\log d\log K)$ 時間となります.
本記事では計算量をオーダーの意味で評価しましたが,算術計算量,すなわち係数環 $R$ における環演算の回数をより詳しく見ると,定数倍にも違いがあります.標準的な高速多項式演算を用いる場合,おおまかには
- Fiduccia のアルゴリズム:$3M(d)\log K$ 程度
- Bostan–Mori のアルゴリズム:$2M(d)\log K$ 程度
と評価でき,この意味では Bostan–Mori のアルゴリズムの方が一般に高速です.さらに FFT を用いる場合には,両者とも定数倍改善の余地がありますが,Bostan–Mori のアルゴリズムでは主要項を $\frac{2}{3}M(d)\log K$ 程度まで小さくできることが知られています.
このような定数倍高速化のテクニックや,本記事で扱わなかった Bostan–Mori のアルゴリズムの発展的な話題については今後の講座で扱います.
また,本記事では線形漸化式あるいは母関数の有理式表示が与えられている数列を扱いました.今回の講座の内容に加えて,線形漸化式の復元について学ぶことで,本記事の手法がさらに競技プログラミングの問題に適用しやすくなります.線形漸化式の復元についても,今後の講座で扱う予定です.
8. 参考文献
- Charles M. Fiduccia.An Efficient Formula for Linear Recurrences.SIAM Journal on Computing.Vol. 14.No. 1.pp. 106–112.1985.DOI: 10.1137/0214007.https://doi.org/10.1137/0214007.
- Alin Bostan, Ryuhei Mori.A Simple and Fast Algorithm for Computing the N-th Term of a Linearly Recurrent Sequence.Proceedings of Symposium on Simplicity in Algorithms(SOSA).pp. 118–132.2021.DOI: 10.1137/1.9781611976496.14.https://doi.org/10.1137/1.9781611976496.14.(Preprint).https://arxiv.org/abs/2008.08822.
- Ryuhei Mori.線形漸化的数列のN項目の計算.Qiita.https://qiita.com/ryuhe1/items/da5acbcce4ac1911f47a.(閲覧日: 2026-05-27).
- Misuki.[tutorial] Bostan-Mori, an elegant algorithm to compute k-th term of linear recurrence in O(dlgdlgk).Codeforces Blog.https://codeforces.com/blog/entry/111862.(閲覧日: 2026-05-27).
- OI Wiki.常系数齐次线性递推.https://oi-wiki.org/math/poly/linear-recurrence/.(閲覧日: 2026-05-27).
トップページ:AtCoder Algorithm Lectures
質問・誤植報告・補足情報など:Discord サーバー