矩阵，递推

以下默认小写字母 \(\mathbf{a},\mathbf{b}\) 等表示的向量是一个行向量，\(\mathbf{1}\) 表示全 \(1\) 行向量，\(\operatorname{diag(\mathbf{a})}\) 表示一个大小 \(|\mathbf{a}|\times |\mathbf{a}|\) 的矩阵，满足 \((i,i)=a_i\) 且其余项为 \(0\)。

缘起：模拟赛题被AI爆标了，原因竟然是：

求 \(\mathbf{a}\times A^{t}\times \mathbf{b^T}\) 的值。

其中 \(A=\mathbf{1^T}\mathbf{c}+\operatorname{diag(d_0,d_1\dots d_{n-1})}\)。

这样特殊的矩阵需要用一些技术处理，下面先了解这些技术。

常系数齐次线性递推

目前只涉及求 \(n\) 阶齐次线性递推数列 \(f\) 的第 \(m\) 项。

一般我们需要解决这样的问题：

\[f_m=\sum_{i=m-n}^{m-1}f_ig_{m-i} \]

给定递推系数 \(g_1\sim g_{n}\) 以及初始项 \(f_0\sim f_{n-1}\)。

观察到这个样子很卷积啊：

\[f_m+\sum_{i=m-n}^{m-1}f_i(-g_{m-i})=0 \]

建立 \(G(z)=1-\sum_{i=1}^nz^ig_i\)，或者说设 \(g'_0=1,g'_i=-g_i\)。

设 \(F(z)=\sum_{i=0}^{\infty}z^if_i\)。

我们要求 \([z^m]F(z)\)。

那么就有：

\[F(z)\times G(z)=R(z) \]

其中 \(R(z)=F(z)\times G(z)(\bmod x^n)\)，因为阶数 \(\ge n\) 的项都没了。

因此，做一次长度为 \(n\) 的多项式乘法就可以求出 \(R(z)\) 了，紧接着：

\[F(z)=\frac{R(z)}{G(z)} \]

Bostan-Mori 算法求远处系数

正是求远处项系数，这是一个用递归理解的算法，每次规模减半。

\[\begin{aligned} [z^m]F(z)&=[z^m]\frac{R(z)}{G(z)}\\ &=[z^m]\frac{R(z)G(-z)}{G(z)G(-z)}\\ &=[z^m]\frac{U_{even}(z^2)+zU_{odd}(z^2)}{V(z^2)}\\ \end{aligned} \]

那么按 \(m\) 的奇偶性保留 \(U_{even}\) 或者 \(U_{odd}\)，问题就只剩下全是平方项的多项式了，因此用 \(z^2\) 代替 \(z\)，问题规模就减半了。

最终剩下上下两个常数，直接算逆元就行。

我们每一轮只会有两次多项式乘法，效率非常高。

 	cin>>m>>n;
    vector<int>f(n+1),g(n+1);
    for(int i=1;i<=n;++i){
        cin>>f[i];f[i]=-f[i];
    }
    f[0]=1;
    for(int i=0;i<n;++i)cin>>g[i];
    g=Mul(f,g);
    while(g.size()>1&&(g.back()==0||g.size()>n))g.pop_back();
    while(m){
        //[x^m]g/f
        vector<int>h=f;
        for(int i=1;i<h.size();i+=2)h[i]=p-h[i];
        vector<int>nf=Mul(h,f),ng=Mul(g,h);
        f.clear();g.clear();
        for(int i=0;i<nf.size();i+=2)f.push_back(nf[i]);    
        for(int i=(m&1);i<ng.size();i+=2)g.push_back(ng[i]);
        m>>=1;
    }
    cout<<g[0]*power(f[0],p-2)%p<<"\n";

另一个办法

我们知道设 \(H(z)\)。

\[H(z)=z^n-g_{1}z^{n-1}-g_{2}z^{n-2}\dots \]

而也就说明：

\[z^n\equiv \sum_{i=0}^{n-1}g_{n-i}z^{i}(\bmod H(z)) \]

这就告诉我们，任意 \(z^{t}\)，无论 \(t\) 有多大，只要 \(t>n\)，我们就能提出一个 \(z^n\) 用后面最高 \(n-1\) 阶的式子代替掉。

则任意 \(z^t\) 都能在模 \(H(z)\) 意义下表达为多项式 \(B(z)\)，其中 \(B(z)\) 阶数不超过 \(n-1\)。

那么假设我们要算 \([z^m]\) 这个系数，可以有：

\[z^m\equiv B(z)(\bmod H(z)) \]

求 \(B(z)\) 可以通过多项式快速幂：就跟普通快速幂差不多，但是将乘法和取模替换为多项式运算。

那么就有：

\[z^m\equiv \sum_{i=0}^{n-1}b_iz^i(\bmod H(z)) \]

那么相应的就有：

\[f_m=\sum_{i=0}^{n-1}f_ib_i \]

同样的可以做到 \(O(n\log n\log m)\)，但是常数非常大。

其本质是 \(a^m\) 不断按递推式拆分拆出的最本源的 \([0,n-1]\) 项的系数。

这个东西有个好处比上面的BM算法优秀，它可以做到 \(f\) 并非标量的计算情况，在下面问题的移步中有阐述。

矩阵快速幂与常系数齐次线性递推的转化

这是一个很重要的观察：若存在多项式 \(F(z)\)，满足 \(F(A)=0\)，那么我们称多项式 \(F\) 为矩阵 \(A\) 的零化多项式。

对于任意的 \(n \times n\) 矩阵 \(A\)，序列 \(S_k = \mathbf{u} A^k \mathbf{v}\) （其中 \(\mathbf{u}\) 是行向量， \(\mathbf{v}\) 是列向量）一定满足一个阶数不超过 \(n\) 的线性递推关系。

这是线性代数中的一个经典结论，主要基于 Cayley-Hamilton 定理。

以下是详细的数学原理和如何找到这个递推。

1. 数学原理：Cayley-Hamilton 定理

定理内容：

设 \(A\) 的特征多项式为 \(p(\lambda) = \det(\lambda I - A)\)。

这是一个 \(n\) 次多项式：\(p(\lambda) = \lambda^n + c_{n-1}\lambda^{n-1} + \dots + c_1\lambda + c_0\)。

Cayley-Hamilton 定理指出，矩阵 \(A\) 满足其自身的特征方程：

\[p(A) = A^n + c_{n-1}A^{n-1} + \dots + c_1 A + c_0 I = 0 \]

推导递推关系：

将上式移项，我们可以表示 \(A^n\)：

\[A^n = - (c_{n-1}A^{n-1} + \dots + c_1 A + c_0 I) \]

如果在等式两边同时乘以 \(A^{k-n}\)（假设 \(k \ge n\)），我们得到：

\[A^k = - \sum_{i=0}^{n-1} c_i A^{k-n+i} \]

这是一个关于矩阵幂次的线性递推。

应用到标量序列：

我们在等式两边左乘向量 \(\mathbf{a}\)，右乘向量 \(\mathbf{1}\)（或者任意向量 \(\mathbf{b}\)）：

\[\mathbf{a} A^k \mathbf{1} = - \sum_{i=0}^{n-1} c_i (\mathbf{a} A^{k-n+i} \mathbf{1}) \]

令 \(S_k = \mathbf{a} A^k \mathbf{1}\)，则有：

\[S_k = - c_{n-1}S_{k-1} - c_{n-2}S_{k-2} - \dots - c_0 S_{k-n} \]

结论：

序列 \(S_k\) 满足常系数线性递推关系，且递推式的特征多项式即为矩阵 \(A\) 的特征多项式（或者是其因式，即最小多项式）。递推的阶数不超过 \(n\)。

关于零化多项式的更多细节，[请移步](零化多项式 - spdarkle - 博客园)