算法 - 齐次线性递推求解和优化及 BM 算法
- \(k\) 项线性递推求第 \(n\) 项的 \(O(k^2\log n)\) 的算法。
- 已知线性递推数列的前若干项，求递推式。
- 结合以上两者可以由线性递推数列的起始若干项得到任意一项。
矩阵快速幂
已知 \(x_n=a_1x_{n-1}+a_2x_{n-2}+\dots+a_kx_{n-k}\) 以及 \(x_0,x_1,\dots,x_{k-1}\)，求 \(x_n\)。
\[ \begin{eqnarray} \left( \begin{array}{c} a_1 & a_2 & a_3 & \ldots & a_{k-1} & a_k \\ 1 & 0 & 0 & \ldots & 0 & 0\\ 0 & 1 & 0 & \ldots & 0 & 0\\ 0 & 0 & 1 & \ldots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & \ldots & 1 & 0 \end{array} \right) \cdot \left( \begin{array}{c} x_{n-1} \\ x_{n-2} \\ x_{n-3} \\ x_{n-4} \\ \vdots \\ x_{n-k} \end{array} \right) = \left( \begin{array}{c} x_{n} \\ x_{n-1} \\ x_{n-2} \\ x_{n-3} \\ \vdots \\ x_{n-k+1} \end{array} \right) \end{eqnarray} \]
剩下的就是矩阵快速幂了。
特征多项式
这部分可略过
已知 \(x_n=a_1x_{n-1}+a_2x_{n-2}+\dots+a_kx_{n-k}\) 以及 \(x_0,x_1,\dots,x_{k-1}\)，求 \(x_n\)。
接下来将通过某种方法将 \(x_{n}\) 的值表示为 \(\alpha_1x_{k-1}+\alpha_2x_{k-2}+\dots+\alpha_kx_0\)，记为 \(\vec x_n=(\alpha_1, \alpha_2,\dots,\alpha_k)\)。
首先显然有 \[ \begin{eqnarray} \vec x_{k}&=&(a_1,a_2,\dots,a_{k-1},a_k) \\ \vec x_{k-1}&=&(1,0,\dots,0,0) \\ \vec x_{k-2}&=&(0,1,\dots,0,0) \\ &\vdots& \\ \vec x_1&=&(0,0,\dots,1,0) \\ \vec x_0&=&(0,0,\dots,0,1) \\ \end{eqnarray} \] 从 \(k+2\) 项开始，利用递推式能够将其表示为前 \(k\) 个向量的线性组合 (\(\vec x_{n+1}=a_1\vec x_{n-1}+a_2\vec x_{n-2}+\dots+a_k\vec x_{n-k}\))。如果得到了 \(x_n\) 的这种表示，那么求出它的值也就轻而易举了。接下来要做的就是加速这个过程。
正文
记 \(x_n\) 的值为 \(\lambda^n\)，问题转化为利用 \({x_n}\) 的特征方程 \(f(\lambda)=-\lambda^k+a_1\lambda^{k-1}+a_2\lambda^{k-2}+\dots+a_k\lambda^0=0\) 以及 \(\lambda^i(0\le i\lt k)\) 的值求 \(\lambda^n\)。
若将 \(\lambda^n\) 表示为 \(g(\lambda)\times f(\lambda)+r(\lambda)\)，其中 \(g(\lambda)\) 是关于 \(\lambda\) 的多项式，\(r(\lambda)\) 是关于 \(\lambda\) 的一个最高次数小于 \(k\) 的多项式。由于 \(f(\lambda)=0\)，所以 \(\lambda^n=r(\lambda)\)，而 \(r(\lambda)\) 的计算就很容易了。
如何求 \(r(\lambda)\) 呢？只需借助多项式快速幂求 \(\lambda^n \mod f(\lambda)\) 即可。
复杂度 \(O(k^2\log n)\)，用 FFT 或 NTT+CRT 可以优化到 \(O(k\log k \log n)\)。
补充
代码量差不多，但效率更差的矩阵快速幂的方法就没有用了吗？如果是双递推（\(\{a_n\}\) 的递推式中有 \(\{b_n\}\)， \(\{b_n\}\) 的递推式中有 \(\{a_n\}\) ）之类的，可能还是要借助矩阵来做。
Berlekamp-Massey (BM) 算法求递推式
已知 \(x_0,x_1,\dots,x_{m}\)，求递推式 \(x_n=a_1x_{n-1}+a_2x_{n-2}+\dots+a_kx_{n-k}\) 在满足数列的同时，\(k\) 尽可能小。
总体思路就是依次考虑每一项，如果发现与当前递推式矛盾则进行修正，同时保证递推项数 \(k\) 尽可能小。
记 \(f(n)=-\lambda^n+a_1\lambda^{n-1}+a_2\lambda^{n-2}+\dots+a_k\lambda^{n-k}(n\ge k)\) （含义与上同），目标就是使 \(\forall n\in[k,m],f(n) =0\) 。
当前值与递推式 \(f(n)\) 矛盾意味着 \(f(u)=d \neq 0\)，若能找到 \(g(n)\) 使 \(g(u)=-1\) 且 \(\forall v\lt u,g(v)=0\)，那么将 \(f(n)\) 修正为 \(f(n) + d\cdot g(n)\) 即可。为了使 \(k\) 尽可能小，\(g\) 的项数要么不超过 \(f\) 的项数，要么尽可能小。
\(g(n)\) 来自于之前的 \(f(n)\) 。对于一次矛盾，得到了 \(f(v)=d\neq 0\) 那么满足 \(g(u)=-1\) 且 \(\forall v\lt u,g(v)=0\) 的一个函数为 \(g(n)=-\frac{f(n)}{d \lambda^{u-v}}\)。
\(g(n)\) 可以重复利用，只是随着 \(u\) 的增长，\(g(n)\) 的项数越来越多。当 \(f(n)\) 的项数小于 \(g(n)\) 时，利用 \(f(n)\) 将 \(g(n)\) 替换掉，这样可以保证 \(k\) 尽可能小。
一开始没有 \(g(n)\) 怎么办呢，将 \(g(n)\) 设成 \(g(n)=C\lambda^n\) 即可，\(C\) 可以是任意数字。
写在最后
代码见模板库。
看了很多网上的资料，几乎都是晦涩难懂。特征多项式这部分我没有涉及到矩阵，应该更好理解了。BM 算法这部分算是对 wiki 的一些补充，但记号可能不太一样。最后还是希望我写的这些对你有所帮助。