矩阵，递推

以下默认小写字母 \(\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\) 的零化多项式。

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

事实上变成了单独开篇

这是一个将一般矩阵快速幂 \(O(n^3\log V)\) 降低为 \(O(n^3+n\log n\log V)\) 的方法。

跟上个链接里的向量情形一致，你只需要求出前 \(n\) 个答案（标量向量都行），然后放 BM 里跑出系数（向量需要随即投影），接着用多项式取模算出系数，再合并答案。