[Tutorial] Bostan-Mori algorithm; compute k-th term of rational power series and linear recurrences

#	User	Rating
1	Benq	3792
2	VivaciousAubergine	3647
3	Kevin114514	3611
4	jiangly	3583
5	strapple	3515
6	tourist	3470
7	Radewoosh	3415
8	Um_nik	3376
9	maroonrk	3361
10	XVIII	3345

#	User	Contrib.
1	Qingyu	162
2	adamant	148
3	Um_nik	146
4	Dominater069	143
5	errorgorn	141
6	cry	138
7	Proof_by_QED	136
8	YuukiS	135
9	chromate00	134
10	soullless	133

Prerequisites

FFT / NTT
Basic understanding of formal power series / generating functions
Basic understanding of linear recurrences

Bostan-Mori algorithm

This algorithm allows us to compute the form $$$[x^k] \dfrac {P(x)} {Q(x)}$$$ where $$$P(x)$$$ and $$$Q(x)$$$ are polynomials with $$$Q(0) \ne 0$$$. It runs in $$$O(d \log d \log k)$$$ time, where $$$d := \max (\deg P, \deg Q)$$$.

The idea for Bostan-Mori is this: we consider a transform $$$\Xi _i: i \in \{ 0, 1 \}$$$ that extracts and reindexes the even-order or odd-order terms of a generating function:

$$$\displaystyle \Xi_i \left(\sum_{n \geq 0} a_n x^n \right) = \sum_{m \geq 0} a_{2m + i} x^m$$$

If we could obtain $$$\Xi_i \left( \dfrac {P(x)} {Q(x)} \right)$$$ from $$$\dfrac {P(x)} {Q(x)} $$$ in time $$$t$$$, we could use the following relation:

$$$\displaystyle \lbrack x^k \rbrack \dfrac {P(x)} {Q(x)} = \lbrack x^{\lfloor k/2 \rfloor} \rbrack \Xi_{k \text{ mod } 2} \left( \dfrac {P(x)} {Q(x)} \right)$$$

to repeatedly halve the order of the desired term, and obtain it in time $$$O(t \log k)$$$.

$$$\Xi_i \left( \dfrac {P(x)} {Q(x)} \right)$$$ can be computed as follows. First multiply $$$Q(-x)$$$ to the numerator and denominator:

$$$\displaystyle \frac{P(x)}{Q(x)} = \frac{P(x)Q(-x)}{Q(x)Q(-x)}$$$

Since $$$Q(x)Q(-x)$$$ is an even polynomial, we may define a polynomial $$$V(x)$$$ such that $$$V(x^2) = Q(x)Q(-x)$$$.

Using this $$$V$$$, we may conclude:

$$$\begin{aligned} \Xi_i\left( \frac{P(x)}{Q(x)} \right) &= \Xi_i\left( \frac{P(x)Q(-x)}{V(x^2)} \right) \\ &= \frac{\Xi_i (P(x) Q(-x))}{V(x)} \end{aligned}$$$

We've thus obtained $$$\Xi_i\left( \dfrac{P(x)}{Q(x)} \right)$$$ in $$$O(d \log d)$$$ time.

Therefore, $$$[x^k] \dfrac {P(x)} {Q(x)}$$$ can be obtained by the following algorithm in $$$O(d \log d \log k)$$$ time:

Input: $$$P(x), Q(x), k$$$
Output: $$$[x^k] \dfrac {P(x)} {Q(x)}$$$
$$$\texttt {while } k \gt 0:$$$
$$$\qquad P(x) \leftarrow \Xi_{k \text{ mod } 2 } (P(x)Q(-x))$$$
$$$\qquad Q(x) \leftarrow \Xi_{0 } (Q(x)Q(-x))$$$
$$$\qquad k \leftarrow \lfloor \frac k 2 \rfloor$$$
$$$\texttt {end while}$$$
$$$\texttt {return } \dfrac {P(0)} {Q(0)}$$$

Problem: Find $$$k$$$-th term of linear recurrence

Problem link: https://judge.yosupo.jp/problem/kth_term_of_linearly_recurrent_sequence

Summary: There is an infinite integer sequence $$$(a_i)$$$ which is a linear recurrence of degree $$$d$$$. We are given its initial terms $$$(a_0, a_1, ..., a_{d-1})$$$, as well as the recurrence coefficients $$$(c_1, c_2, ..., c_d)$$$, such that for $$$i \ge d$$$, $$$a_i = \sum_{j=1}^{d} c_j a_{i-j}$$$. For some integer $$$k \leq 10^{18}$$$, find $$$a_k$$$.

Solution

Extra problem: Fibonacci Revisited

Problem link: https://atcoder.jp/contests/abc300/tasks/abc300_h

Summary: We have a linear recurrence where $$$(a_0, a_1, ..., a_{d-1})$$$ and $$$(c_1, c_2, ..., c_d)$$$ are all $$$1$$$, but now we need to find $$$\displaystyle \sum_{m \text{ & } k = k} [x^m]A(x)$$$ instead, where $$$\text{&}$$$ denotes the bitwise AND operation.

Solution