Hi codeforces! I want to show that how to implement the algorithm from [user:noshi91,2024-03-29] who showed that we could modify the Bostan--Mori's algorithm to make use in FPS composition and compositional inverse.↵
↵
[user:noshi91,2024-03-29]'s blog: [FPS の合成と逆関数、冪乗の係数列挙 Θ(n (log(n))^2)](https://noshi91.hatenablog.com/entry/2024/03/16/224034) in Japanese.↵
↵
## Compositional Inverse↵
↵
The key idea is consider the bivariate formal power series↵
↵
$$↵
\dfrac{1}{1-yf(x)}\in\mathbb{C}\left\lbrack\left\lbrack x,y\right\rbrack\right\rbrack↵
$$↵
↵
Let $f(0)=0$, we have↵
↵
$$↵
\begin{aligned}↵
\left\lbrack x^k\right\rbrack\frac{1}{1-yf(x)}=\left\lbrack x^k\right\rbrack\sum_{i\geq 0}y^if(x)^i↵
\end{aligned}↵
$$↵
↵
Apply Bostan--Mori's algorithm here:↵
↵
$$↵
\begin{aligned}↵
\frac{P(x,y)}{Q(x,y)}&=\frac{P(x,y)Q(-x,y)}{Q(x,y)Q(-x,y)} \\↵
&=\frac{U_{\text{even}}(x^2,y)+xU_{\text{odd}}(x^2,y)}{V(x^2,y)}↵
\end{aligned}↵
$$↵
↵
so↵
↵
$$↵
\left\lbrack x^k\right\rbrack\frac{P(x,y)}{Q(x,y)}=↵
\begin{cases}↵
\left\lbrack x^{k/2}\right\rbrack &\dfrac{U_{\text{even}}(x,y)}{V(x,y)},\space &\text{if }k\equiv 0\pmod{2} \\↵
\left\lbrack x^{(k-1)/2}\right\rbrack &\dfrac{U_{\text{odd}}(x,y)}{V(x,y)},\space &\text{if }k\equiv 1\pmod{2}↵
\end{cases}↵
$$↵
↵
The degree of $y$ doubled after each time of iteration, but the degree of $x$ halved.↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Enum-}k\operatorname{th-Term-of-Power}(f,k,n)\text{:} \\↵
&\textbf{Input}\text{: } f\in x\mathbb{C}\left\lbrack x\right\rbrack,k,n\in\mathbb{Z}_{\geq 0}\text{.} \\↵
&\textbf{Output}\text{: }\left(\left\lbrack x^k\right\rbrack f^0,\left\lbrack x^k\right\rbrack f,\left\lbrack x^k\right\rbrack f^2,\dots ,\left\lbrack x^k\right\rbrack f^{n-1}\right)\text{.} \\↵
1&P(x,y)\gets 1 \\↵
2&Q(x,y)\gets 1-yf \\↵
3&\textbf{while }\deg_y(P)<n-1\textbf{ or }k\neq 0\textbf{ do}\text{:} \\↵
4&\qquad U_0(x^2,y)+xU_1(x^2,y)\gets P(x,y)Q(-x,y) \\↵
5&\qquad V(x^2,y)\gets Q(x,y)Q(-x,y) \\↵
6&\qquad P\gets U_{k\bmod{2}}(x,y) \pbmod{x^{\left\lfloor k/2\right\rfloor +1}} \\↵
7&\qquad Q\gets V(x,y) \pbmod{x^{\left\lfloor k/2\right\rfloor +1}} \\↵
8&\qquad k\gets \left\lfloor k/2\right\rfloor \\↵
9&\textbf{end while} \\↵
10&A(y)=\sum_{j=0}^{\deg_y(P)}a_jy^j\gets \left\lbrack x^0\right\rbrack P(x,y) \\↵
11&\textbf{return }(a_0,\dots ,a_{n-1})↵
\end{array}↵
$$↵
↵
Combined with Lagrange Inversion Formula, we have the following algorithm:↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Compositional-Inverse}(f,n)\text{:} \\↵
&\textbf{Input}\text{: }f\in x\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack, f'(0)\neq 0,n\in\mathbb{Z}_{\geq 2}\text{.} \\↵
&\textbf{Output}\text{: }g\pmod{x^n} \text{ such that }f(g)\equiv g(f)\equiv x\pmod{x^n}\text{.} \\↵
1&t\gets f'(0) \\↵
2&F(x)\gets f\left(t^{-1}x\right) \\↵
3&(a_0,\dots ,a_{n-1})\gets \operatorname{Enum-}k\operatorname{th-Term-of-Power}(F,n-1,n) \\↵
&\because \left\lbrack x^{n-1}\right\rbrack F^k=\frac{k}{n-1}\left\lbrack x^{n-1-k}\right\rbrack \left(F^{\langle -1\rangle}/x\right)^{-(n-1)} \\↵
&\therefore \frac{n-1}{k}a_{k}=\left\lbrack x^{n-1-k}\right\rbrack\left(F^{\langle -1\rangle}/x\right)^{-(n-1)},\space (k>0) \\↵
4&G(x)\gets \sum_{k=1}^{n-1}\frac{n-1}{k}a_{k}x^{n-1-k} \\↵
5&F^{\langle -1\rangle}/x\gets \left(G(x)^{1/(n-1)}\right)^{-1}\mod{x^{n-1}} \\↵
&\because (F^{\langle -1\rangle}/x)(0)=1 \\↵
&\therefore F^{\langle -1\rangle}/x=\exp\left(\frac{1}{1-n}\log G\right) \\↵
&\because F^{\langle -1\rangle}\circ f\circ \left(t^{-1}x\right)=x \\↵
&\therefore \left((t^{-1}x) \circ F^{\langle -1\rangle}\right)\circ f=\left((t^{-1}x) \circ F^{\langle -1\rangle}\right)\circ f \circ \left((t^{-1}x)\circ (tx)\right)=x \\↵
6&\textbf{return }\left((t^{-1}x) \circ F^{\langle -1\rangle}\right)↵
\end{array}↵
$$↵
↵
## Composition↵
↵
Given $f=\sum_{k\geq 0}f_kx^k\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack,g\in x\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack$. then↵
↵
$$↵
f(g)=\sum_{k\geq 0}f_kg^k↵
$$↵
↵
Consider again↵
↵
$$↵
\dfrac{f(y^{-1})}{1-yg(x)}=\sum_{k\geq 0}(\cdots +f_ky^{-k}+\cdots)y^kg(x)^k↵
$$↵
↵
Our goal is computing↵
↵
$$↵
\left\lbrack y^0\right\rbrack \dfrac{f(y^{-1})}{1-yg(x)}↵
$$↵
↵
Apply Bostan--Mori's algorithm here:↵
↵
$$↵
\begin{aligned}↵
\frac{P(y)}{Q(x,y)}&=\frac{P(y)}{Q(x,y)Q(-x,y)}Q(-x,y) \\↵
&=\frac{P(y)}{V(x^2,y)}Q(-x,y)↵
\end{aligned}↵
$$↵
↵
We want to compute $\dfrac{P(y)}{Q(x,y)Q(-x,y)}$ first, and we only need to track finitely many terms of $y^j$ for $j\leq 0$.↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Composition-Subprocedure}(P,Q,n)\text{:} \\↵
&\textbf{Input}\text{: }P=\sum_{0\leq j\leq n}p_jy^{-j}\in\mathbb{C}((y)),Q\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack\left\lbrack y\right\rbrack ,Q(0)=1,[x^{\lbrack 0,k)}y^k]Q=(0,\dots ,0),(\forall k>0)\text{.} \\↵
&\textbf{Output}\text{: }\left\lbrack y^{\left\lbrack -\deg_y(Q)+1,0\right\rbrack}\right\rbrack\dfrac{P}{Q}\bmod{x^{n+1}}\text{.} \\↵
1&d\gets \deg_y(Q)\\↵
2&\textbf{if }n=0\textbf{ then return }\left(\left\lbrack y^{-d+1}\right\rbrack P,\dots ,\left\lbrack y^0\right\rbrack P\right) \\↵
&n=0\iffmplies \deg_x(Q)=0 \\↵
&\because Q(0)=1,[x^{\lbrack 0,k)}y^k]Q=(0,\dots ,0),(\forall k>0) \\↵
&\therefore Q\equiv 1\pmod{x^1},\space \left(\left\lbrack y^{-d+1}\right\rbrack P,\dots ,\left\lbrack y^0\right\rbrack P\right)=\left(\left\lbrack y^{-d+1}\right\rbrack P/Q,\dots ,\left\lbrack y^0\right\rbrack P/Q\right) \\↵
3&V(x^2,y)\gets Q(x,y)Q(-x,y)\bmod{x^{n+1}} \\↵
4&(t_{-2d+1},\dots ,t_0)\gets \operatorname{Composition-Subprocedure}\left(P,V(x,y),\left\lfloor n/2\right\rfloor\right) \\↵
5&T(x,y)\gets \sum_{j=-2d+1}^0t_jy^j \\↵
6&U(x,y)=\sum_{j=-2d+1}^d u_jy^j\gets T(x^2,y)Q(-x,y)\bmod{x^{n+1}} \\↵
7&\textbf{return }\left(u_{-d+1},\dots ,u_0\right)↵
\end{array}↵
$$↵
↵
and the composition algorithm is just↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Composition}(f,g,n)\text{:} \\↵
&\textbf{Input}\text{: }f\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack,g\in\mathbb{C}\left\lbrack x\right\rbrack ,n\in\mathbb{Z}_{>0}\text{.} \\↵
&\textbf{Output}\text{: }f(g)\bmod{x^n}\text{.} \\↵
1&c\gets g(0) \\↵
2&F(x)\gets f(x+c)\bmod{x^n} \\↵
3&G(x)\gets (g-c)\bmod{x^n} \\↵
4&(u)\gets\operatorname{Composition-Subprocedure}\left(F\left(y^{-1}\right),1-yG,n-1\right) \\↵
&f(g)=f\circ (x+c)\circ G=F(G) \\↵
5&\textbf{return }(u)↵
\end{array}↵
$$↵
↵
One can tryonwith [composition_of_formal_power_series_large](https://judge.yosupo.jp/problem/composition_of_formal_power_series_large) at [Library Checker](https://judge.yosupo.jp/).↵
↵
## Reference↵
- Alin Bostan, Ryuhei Mori. [A Simple and Fast Algorithm for Computing the N-th Term of a Linearly Recurrent Sequence](https://arxiv.org/abs/2008.08822).↵
- [user:noshi91,2024-03-29]. [FPS の合成と逆関数、冪乗の係数列挙 Θ(n (log(n))^2)](https://noshi91.hatenablog.com/entry/2024/03/16/224034).
↵
[user:noshi91,2024-03-29]'s blog: [FPS の合成と逆関数、冪乗の係数列挙 Θ(n (log(n))^2)](https://noshi91.hatenablog.com/entry/2024/03/16/224034) in Japanese.↵
↵
## Compositional Inverse↵
↵
The key idea is consider the bivariate formal power series↵
↵
$$↵
\dfrac{1}{1-yf(x)}\in\mathbb{C}\left\lbrack\left\lbrack x,y\right\rbrack\right\rbrack↵
$$↵
↵
Let $f(0)=0$, we have↵
↵
$$↵
\begin{aligned}↵
\left\lbrack x^k\right\rbrack\frac{1}{1-yf(x)}=\left\lbrack x^k\right\rbrack\sum_{i\geq 0}y^if(x)^i↵
\end{aligned}↵
$$↵
↵
Apply Bostan--Mori's algorithm here:↵
↵
$$↵
\begin{aligned}↵
\frac{P(x,y)}{Q(x,y)}&=\frac{P(x,y)Q(-x,y)}{Q(x,y)Q(-x,y)} \\↵
&=\frac{U_{\text{even}}(x^2,y)+xU_{\text{odd}}(x^2,y)}{V(x^2,y)}↵
\end{aligned}↵
$$↵
↵
so↵
↵
$$↵
\left\lbrack x^k\right\rbrack\frac{P(x,y)}{Q(x,y)}=↵
\begin{cases}↵
\left\lbrack x^{k/2}\right\rbrack &\dfrac{U_{\text{even}}(x,y)}{V(x,y)},\space &\text{if }k\equiv 0\pmod{2} \\↵
\left\lbrack x^{(k-1)/2}\right\rbrack &\dfrac{U_{\text{odd}}(x,y)}{V(x,y)},\space &\text{if }k\equiv 1\pmod{2}↵
\end{cases}↵
$$↵
↵
The degree of $y$ doubled after each time of iteration, but the degree of $x$ halved.↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Enum-}k\operatorname{th-Term-of-Power}(f,k,n)\text{:} \\↵
&\textbf{Input}\text{: } f\in x\mathbb{C}\left\lbrack x\right\rbrack,k,n\in\mathbb{Z}_{\geq 0}\text{.} \\↵
&\textbf{Output}\text{: }\left(\left\lbrack x^k\right\rbrack f^0,\left\lbrack x^k\right\rbrack f,\left\lbrack x^k\right\rbrack f^2,\dots ,\left\lbrack x^k\right\rbrack f^{n-1}\right)\text{.} \\↵
1&P(x,y)\gets 1 \\↵
2&Q(x,y)\gets 1-yf \\↵
3&\textbf{while }\deg_y(P)<n-1\textbf{ or }k\neq 0\textbf{ do}\text{:} \\↵
4&\qquad U_0(x^2,y)+xU_1(x^2,y)\gets P(x,y)Q(-x,y) \\↵
5&\qquad V(x^2,y)\gets Q(x,y)Q(-x,y) \\↵
6&\qquad P\gets U_{k\bmod{2}}(x,y) \
7&\qquad Q\gets V(x,y) \
8&\qquad k\gets \left\lfloor k/2\right\rfloor \\↵
9&\textbf{end while} \\↵
10&A(y)=\sum_{j=0}^{\deg_y(P)}a_jy^j\gets \left\lbrack x^0\right\rbrack P(x,y) \\↵
11&\textbf{return }(a_0,\dots ,a_{n-1})↵
\end{array}↵
$$↵
↵
Combined with Lagrange Inversion Formula, we have the following algorithm:↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Compositional-Inverse}(f,n)\text{:} \\↵
&\textbf{Input}\text{: }f\in x\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack, f'(0)\neq 0,n\in\mathbb{Z}_{\geq 2}\text{.} \\↵
&\textbf{Output}\text{: }g\pmod{x^n} \text{ such that }f(g)\equiv g(f)\equiv x\pmod{x^n}\text{.} \\↵
1&t\gets f'(0) \\↵
2&F(x)\gets f\left(t^{-1}x\right) \\↵
3&(a_0,\dots ,a_{n-1})\gets \operatorname{Enum-}k\operatorname{th-Term-of-Power}(F,n-1,n) \\↵
&\because \left\lbrack x^{n-1}\right\rbrack F^k=\frac{k}{n-1}\left\lbrack x^{n-1-k}\right\rbrack \left(F^{\langle -1\rangle}/x\right)^{-(n-1)} \\↵
&\therefore \frac{n-1}{k}a_{k}=\left\lbrack x^{n-1-k}\right\rbrack\left(F^{\langle -1\rangle}/x\right)^{-(n-1)},\space (k>0) \\↵
4&G(x)\gets \sum_{k=1}^{n-1}\frac{n-1}{k}a_{k}x^{n-1-k} \\↵
5&F^{\langle -1\rangle}/x\gets \left(G(x)^{1/(n-1)}\right)^{-1}\mod{x^{n-1}} \\↵
&\because (F^{\langle -1\rangle}/x)(0)=1 \\↵
&\therefore F^{\langle -1\rangle}/x=\exp\left(\frac{1}{1-n}\log G\right) \\↵
&\because F^{\langle -1\rangle}\circ f\circ \left(t^{-1}x\right)=x \\↵
&\therefore \left((t^{-1}x) \circ F^{\langle -1\rangle}\right)\circ f=\left((t^{-1}x) \circ F^{\langle -1\rangle}\right)\circ f \circ \left((t^{-1}x)\circ (tx)\right)=x \\↵
6&\textbf{return }\left((t^{-1}x) \circ F^{\langle -1\rangle}\right)↵
\end{array}↵
$$↵
↵
## Composition↵
↵
Given $f=\sum_{k\geq 0}f_kx^k\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack,g\in x\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack$
↵
$$↵
f(g)=\sum_{k\geq 0}f_kg^k↵
$$↵
↵
Consider again↵
↵
$$↵
\dfrac{f(y^{-1})}{1-yg(x)}=\sum_{k\geq 0}(\cdots +f_ky^{-k}+\cdots)y^kg(x)^k↵
$$↵
↵
Our goal is computing↵
↵
$$↵
\left\lbrack y^0\right\rbrack \dfrac{f(y^{-1})}{1-yg(x)}↵
$$↵
↵
Apply Bostan--Mori's algorithm here:↵
↵
$$↵
\begin{aligned}↵
\frac{P(y)}{Q(x,y)}&=\frac{P(y)}{Q(x,y)Q(-x,y)}Q(-x,y) \\↵
&=\frac{P(y)}{V(x^2,y)}Q(-x,y)↵
\end{aligned}↵
$$↵
↵
We want to compute $\dfrac{P(y)}{Q(x,y)Q(-x,y)}$ first, and we only need to track finit
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Composition-Subprocedure}(P,Q,n)\text{:} \\↵
&\textbf{Input}\text{: }P=\sum_{0\leq j\leq n}p_jy^{-j}\in\mathbb{C}((y)),Q\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack\left\lbrack y\right\rbrack ,Q(0)=1,[x^{\lbrack 0,k)}y^k]Q=(0,\dots ,0),(\forall k>0)\text{.} \\↵
&\textbf{Output}\text{: }\left\lbrack y^{\left\lbrack -\deg_y(Q)+1,0\right\rbrack}\right\rbrack\dfrac{P}{Q}\bmod{x^{n+1}}\text{.} \\↵
1&d\gets \deg_y(Q)\\↵
2&\textbf{if }n=0\textbf{ then return }\left(\left\lbrack y^{-d+1}\right\rbrack P,\dots ,\left\lbrack y^0\right\rbrack P\right) \\↵
&n=0\i
&\because Q(0)=1,[x^{\lbrack 0,k)}y^k]Q=(0,\dots ,0),(\forall k>0) \\↵
&\therefore Q\equiv 1\pmod{x^1},\space \left(\left\lbrack y^{-d+1}\right\rbrack P,\dots ,\left\lbrack y^0\right\rbrack P\right)=\left(\left\lbrack y^{-d+1}\right\rbrack P/Q,\dots ,\left\lbrack y^0\right\rbrack P/Q\right) \\↵
3&V(x^2,y)\gets Q(x,y)Q(-x,y)\bmod{x^{n+1}} \\↵
4&(t_{-2d+1},\dots ,t_0)\gets \operatorname{Composition-Subprocedure}\left(P,V(x,y),\left\lfloor n/2\right\rfloor\right) \\↵
5&T(x,y)\gets \sum_{j=-2d+1}^0t_jy^j \\↵
6&U(x,y)=\sum_{j=-2d+1}^d u_jy^j\gets T(x^2,y)Q(-x,y)\bmod{x^{n+1}} \\↵
7&\textbf{return }\left(u_{-d+1},\dots ,u_0\right)↵
\end{array}↵
$$↵
↵
and the composition algorithm is just↵
↵
$$↵
\begin{array}{ll}↵
&\textbf{Algorithm }\operatorname{Composition}(f,g,n)\text{:} \\↵
&\textbf{Input}\text{: }f\in\mathbb{C}\left\lbrack\left\lbrack x\right\rbrack\right\rbrack,g\in\mathbb{C}\left\lbrack x\right\rbrack ,n\in\mathbb{Z}_{>0}\text{.} \\↵
&\textbf{Output}\text{: }f(g)\bmod{x^n}\text{.} \\↵
1&c\gets g(0) \\↵
2&F(x)\gets f(x+c)\bmod{x^n} \\↵
3&G(x)\gets (g-c)\bmod{x^n} \\↵
4&(u)\gets\operatorname{Composition-Subprocedure}\left(F\left(y^{-1}\right),1-yG,n-1\right) \\↵
&f(g)=f\circ (x+c)\circ G=F(G) \\↵
5&\textbf{return }(u)↵
\end{array}↵
$$↵
↵
One can try
↵
## Reference↵
- Alin Bostan, Ryuhei Mori. [A Simple and Fast Algorithm for Computing the N-th Term of a Linearly Recurrent Sequence](https://arxiv.org/abs/2008.08822).↵
- [user:noshi91,2024-03-29]. [FPS の合成と逆関数、冪乗の係数列挙 Θ(n (log(n))^2)](https://noshi91.hatenablog.com/entry/2024/03/16/224034).
