线性递推与整式递推数列

线性递推与整式递推数列

本文主要摘录自 2019 年国家候选队论文集《两类递推数列的性质和应用》——钟子谦。

线性递推数列

定义

对于无限数列 \(\{ a_0,a_1,\cdots\}\)，和有限非空数列 \(\{ r_0,r_1,\cdots,r_{m-1} \}\)，若对于任意 \(p\geq m-1\)，有 \(\sum_{k=0}^{m-1}a_{p-k}r_k=0\)，则称 \(r\) 为 \(a\) 的线性递归式。特别的，若 \(r_0=1\)，称 \(r\) 为 \(a\) 的线性递推式。

我们称存在线性递推式的无限数列为线性递推数列。

对于有限数列，我们类似的定义线性递推式，称线性递推式的阶数，为其长度减一（上文中的 \(m-1\) 为 \(r\) 的阶数），称数列 \(a\) 阶数最小的线性递推式为 \(a\) 的最短线性递推式。

基本性质与判定

从生成函数的角度来看线性递归式，有结论：

对于无限数列 \(a\)，和有限非空数列 \(r\)，记两者的生成函数为 \(A,R\)，数列 \(r\) 为 \(a\) 的线性递归式等价于存在次数不超过 \(m-2\) 的多项式 \(S\)，满足 \(AR+S=0\)。

有限数列的情况仅需改成 \(AR+S\equiv 0\pmod {x^n}\)。

结论的证明很显然，考虑其推论：

对于无限数列 \(a\)，记其生成函数为 \(A\)，\(a\) 为线性递推数列等价于存在常数项为 \(1\) 的多项式 \(R\) 和多项式 \(S\) 满足 \(A=\dfrac S R\)。数列 \(a\) 的最短递推式阶数为所有的 \(R,S\) 对中 \(\max(\deg R,\deg S+1)\) 的最小可能值。

线性递推数列满足以下的封闭性：

1. 线性递推数列乘常数仍是线性递推数列。
2. 线性递推数列平移仍是线性递推数列。
3. 两个线性递推数列相加仍是线性递推数列。
4. 两个线性递推数列卷积仍是线性递推数列。
5. 两个线性递推数列点乘仍是线性递推数列。

有关算法

求出数列的最短线性递推式

我们先考虑有限数列的情况。
一种简单的做法是高斯消元解出最短线性递推式，对于阶为 \(n\) 的有限数列来说复杂度 \(O(n^3)\)，使用 Berlekamp-Massey(BM) 算法可以将复杂度降至 \(O(n^2)\)。（实际上其复杂度为 \(O(nm)\)，\(m\) 为最短线性递推式阶数。）

BM 会对数列 \(a\) 的每个前缀求出其最短线性递推式。

为了方便，定义新的递推系数的形式 \(\{ f_0,f_1\cdots,f_{m-1}\}\)：
\(a_i=\sum_{j=0}^{m-1}f_ja_{i-1-j},\forall i\geq m\)。上文定义的阶数在新的形式下即为 \(|f|\)。

按位考虑 \(a\) 的每一位，在递推出现错误时对系数进行调整，记 \(\{a_0,a_1,\cdots,a_i\}\) 的最短线性递推式为 \(F_i\)，记 \(l_i=|F_i|\)。

初始时 \(F_{-1}=\varnothing\)，处理到 \(a_i\) 时有两种情况：

1. \(F_{i-1}\) 对 \(a_i\) 也成立，不需要调整，\(F_i=F_{i-1}\)。
2. \(F_{i-1}\) 对 \(a_i\) 不成立，考虑进行调整。

记 \(\Delta_i=a_i-\sum_{j=0}F_{i-1,j}a_{i-1-j}\)。

若此时是第一次对递推式进行修改（\(a_i\) 为第一个非 0 项），直接令 \(F_i\) 为 \(i+1\) 个 \(0\)。
否则找到上一次对递推式修改的位置 \(k\)，考虑构造一个序列 \(G=\{g_0,g_1,\cdots,g_{m-1}\}\)，满足：\(\sum_{j=0}^{m-1}g_ja_{p-1-j}=0,\forall p\in [m,i-1]\) 且 \(\sum_{j=0}^{m-1}g_ja_{i-1-j}=\Delta_i\)。
那么我们令 \(F_i=F_{i-1}+G\) 即可。

一种可行的构造是 \(G=\{0,0,\cdots,0,\dfrac{\Delta_i}{\Delta_k},-\dfrac{\Delta_i}{\Delta_k}F_{k-1}\}\)。
前缀 0 的个数为 \(i-k-1\)，最后一项的意思是将 \(F_{k-1}\) 乘上系数加在 \(G\) 后。
不难验证构造的正确性，我们还需说明构造的最小性。

引理：若 \(F_{i-1}\) 不是 \(\{a_0,a_1,\cdots,a_i\}\) 的最短线性递推式，有：

\[l_i\geq \max(l_{i-1},i+1-l_{i-1}) \]

考虑反证，假设 \(l_i\leq i-l_{i-1}\)。设 \(F_{i-1}=\{p_j\},F_i=\{q_j\}\)。

\[\sum_{j=0}^{l_{i-1}-1}p_ja_{i-1-j}=\sum_{j=0}^{l_{i-1}-1}p_j\sum_{k=0}^{l_i-1}q_ka_{i-1-j-1-k}\\=\sum_{k=0}^{l_i-1}q_k\sum_{j=0}^{l_{i-1}-1}p_ja_{i-1-k-1-j}=\sum_{k=0}^{l_i-1}q_ka_{i-1-k}=a_i \]

那么 \(F_{i-1}\) 为 \(\{a_0,a_1,\cdots,a_i\}\) 的最短递推式，矛盾。

回到我们上文给出的构造，当 \(F_{i-1}\) 不是\(\{a_0,a_1,\cdots,a_i\}\) 的最短递推式时，我们构造得到的 \(l_i\leftarrow l_{k-1}+i-k\)。

考虑我们归纳这样的 \(l_i\) 就是最短的，那么有 \(l_k=k+1-l_{k-1}\)，则 \(l_{i-1}=k+1-l_{k-1}\)。故新的 \(l_i=i+1-l_{i-1}\)。构造的最小性成立。

对于无限序列，若我们知道其最短递推式的长度 \(m\)，我们取其前 \(2m\) 项跑 BM 就能得到其最短递推式。

实现很清新：

vec BM(int *a,int n){
    int lst=-1,val=0;
    vec F0,F1;F0.clear();F1.clear();
    for(int i=0;i<=n;i++){
        int del=a[i];
        for(int j=0;j<F1.size();j++)del=sub(del-1ll*F1[j]*a[i-j-1]%mod);
        if(del==0)continue;
        if(lst==-1){lst=i,val=del;F0=F1;for(int j=0;j<=i;j++)F1.eb(0);continue;}
        vec G;G.clear();
        for(int j=1;j<=i-lst-1;j++)G.eb(0);
        int w=1ll*del*power(val,mod-2)%mod;
        G.eb(w);for(int j=0;j<F0.size();j++)G.eb(sub(-1ll*w*F0[j]%mod));
        F0=F1;F1.resize(max(F1.size(),G.size()));
        for(int j=0;j<G.size();j++)F1[j]=add(F1[j]+G[j]);
        lst=i,val=del;
    }
    return F1;
}

求出一个线性递推数列的第 \(n\) 项

常系数齐次线性递推，运用 BM 得到线性递推式后，使用多项式取模以 \(O(m^2\log n)\) 或者 \(O(m\log m\log n)\) 的复杂度求得数列第 \(n\) 项。

应用

求向量序列和矩阵序列的最短递推式

考虑如何求出 \(n\) 维行向量序列 \(\{t_0,t_1,\cdots\}\) 的线性递推式。假设考虑在模 \(p\) 意义下随机一个 \(n\) 维列向量 \(v\)，转而计算 \(\{t_0v,t_1v,\cdots\}\) 这个标量序列的最短线性递推式。
根据著名的 Schwartz-Zippel 引理，我们可以推导出至少有 \(1-\dfrac n p\) 的概率我们得到的线性递推式正确。
我们类似的求 \(n\times m\) 的矩阵序列 \(\{t_0,t_1,\cdots\}\) 的线性递推式，随机一个 \(n\) 维列向量 \(u\)，和一个 \(m\) 维行向量，计算 \(\{ut_0v,ut_1v,\cdots\}\) 的线性递推式。我们有至少 \(1-\dfrac{n+m}p\) 的概率得到正确线性递推式。

求矩阵的最小多项式

定义 \(n\times n\) 的矩阵 \(M\) 的最小多项式为次数最小的使得 \(f(M)=0\) 的多项式 \(f\)。
矩阵的最小多项式显然就是 \(\{I,M,M^2,\cdots\}\) 的最短线性递推式，我们使用上文的方法求出矩阵序列的最短线性递推式即可。
具体的，我们可以在 \(O(n^3)\) 的时间内将标量序列得到，对于稀疏矩阵，记 \(e\) 为矩阵中的非零元素个数，可以在 \(O(n(n+e))\) 的时间内得到标量序列。
之后进行最劣 \(O(n^2)\) 的 BM 即可。

优化动态规划

传统使用矩阵快速幂优化的 \(F_n=F_0M^n\) 形式的动态规划复杂度为 \(O(m^3\log n)\)。
由于 \(\{F_i\}\) 为线性递推向量序列，我们使用 BM 求出其线性递推式，再使用常系数齐次线性递推，可以做到 \(O(m^3+m^2\log n)\) 或者 \(O(m^3+m\log m\log n)\)。

解稀疏线性方程组

即求 \(x=A^{-1}b\)。
考虑求得 \(\{b,Ab,A^2b,\cdots\}\) 的最短递推式 \(r\)，我们有 \(\sum_{i=0}^{m-1}A^ibr_{m-1-i}=0\)，由于 \(r\) 最短，\(r_{m-1}\neq 0\)。在两边乘上 \(A^{-1}\) 移项，可以得到 \(A^{-1}b=-\dfrac 1 {r_{m-1}}(\sum_{i=0}^{m-2}A^ibr_{m-2-i})\)。
瓶颈在于求 \(\{b,Ab,A^2b,\cdots,A^{2n}b\}\)，若 \(A\) 中有 \(e\) 个非零元素，总时间复杂度 \(O(n(n+e))\)。

求稀疏矩阵行列式

注意到我们可以快速地求出稀疏矩阵的最小多项式，而当矩阵的每个特征值的几何重数均为一时最小多项式就是特征多项式。对于 \(n\) 阶矩阵，特征多项式的常数项乘上 \((-1)^n\) 即为行列式（因为行列式即全部特征值的乘积）。

由于我们所求的矩阵 \(A\) 不一定每个特征值的几何重数均为一，我们给它乘上一个随机的对角矩阵 \(B\)，可以证明 \(AB\) 有至少 \(1-\dfrac{2n^2-n} p\) 的概率满足性质，得到 \(\det(AB)\) 后由于 \(\det(AB)=\det A\det B\)，\(\det B\) 容易计算，我们容易得到 \(\det A\) 的值。

求稀疏矩阵的秩

感觉没啥用啊，有点复杂，鸽着。

整式递推数列

介于笔者的高数水平一般，这个部分很多会省略证明过程。

定义

对于无限数列 \(a\) 和有限非空多项式列 \(\{P_0,P_1,\cdots,P_{m-1}\}\)，若 \(P_0\) 非 \(0\) 且对于任意 \(p\geq m-1\)，有 \(\sum_{k=0}^{m-1}a_{p-k}P_k(p)=0\)，则称 \(P\) 为数列 \(a\) 的整式递推式，称存在整式递推式的数列 \(a\) 为整式递推数列。

对于有限数列的定义类似。

显然，上述的线性递推数列是整式递推数列。
我们类比线性递推，定义整式递推式的阶数为 \(m-1\)，次数为 \(\max_{i=0}^{m-1}\deg P_i\)。

为了更好的描述整式递推数列的性质、连接生成函数，我们引入微分有限（D-finite）和代数形式幂级数（Algebraic）。
我们称一个形式幂级数 \(A(x)\) 为微分有限的当且仅当存在多项式数列 \(\{Q_0(x),Q_1(x),\cdots,Q_{m-1}(x)\}\)，满足 \(Q_{m-1}\neq 0\) 且 \(\sum_{i=0}^{m-1}Q_i(x)A^{(i)}(x)=0\)。
我们称一个形式幂级数 \(A(x)\) 为代数形式幂级数当且仅当其在 \(\mathbb K(x)\) 上是代数的（即 \(A(x)\) 为系数在 \(\mathbb K(x)\) 内的多项式方程的根）。

代数形式幂级数的定义类似于数域中代数数的定义，与其相对的为超越数。
举一些容易理解的例子，代数数：\(12,\sqrt 5\)，超越数：\(\pi,e,\sin 1\)。

基本性质与判定

记一个数列 \(\{a_0,a_1,a_2,\cdots\}\) 的 OGF 为 \(A(x)\)，那么 \(a\) 是整式递推数列当且仅当 \(A\) 微分有限。

下文我们考虑一种常见的生成函数——代数形式幂级数。我们可以证明代数形式幂级数都是微分有限的，不仅如此，我们有更强的结论：

代数形式幂级数 \(u\) 的微分 \(u'\) 也是代数形式幂级数。

记 \(f\) 微分有限，\(g\) 是代数的，有定义的（\(f\) 有限或 \(g\) 常数项为 \(0\)） \(f\circ g\) 微分有限。

与线性递推类似，整式递推数列也满足以下封闭性：

1. 整式递推数列乘常数仍为整式递推数列。
2. 整式递推数列平移仍为整式递推数列。
3. 两个整式递推数列相加仍为整式递推数列。
4. 两个整式递推数列卷积仍为整式递推数列。
5. 两个整式递推数列点乘仍为整式递推数列。

有关算法

求出数列的整式递推式

一个暴力但实用的方法是将 \(\text{阶数}\times \text{次数}\) 个系数直接高斯消元，这里的阶数和次数都是一个预估值，取肯定比真正阶数与次数略大的值即可。

求出一个整式递推数列的第 \(n\) 项

记整式递推式的次数为 \(d\)，阶数为 \(m\)，那么我们直接递推，有暴力的 \(O(nmd)\) 做法。

对于 \(n\) 远大于 \(m,d\) 的情形，存在一个 \(O(\sqrt{nd}(m^3+m^2\log(nd)))\) 的做法。
当 \(m\) 为常数时，其复杂度为 \(O(\sqrt{nd}\log(nd))\)。