矩阵快速幂模板题题解

给你一个矩阵，求矩阵幂，这是矩阵快速幂模板题。但是指数特别大，做不了，G 了。咋办？

引入 Cayley–Hamilton 定理：对于矩阵 \(A\) 的特征多项式 \(p_A(\lambda) = |\lambda E - A|\)，有 \(p_A(A) = \mathbf{0}\)。想看证明的可以翻到题解最底下。

那这个定理对这个题有什么作用呢？容易发现 \(A\) 的特征多项式是 \(k\) 次的，但是 \(A^n\) 是 \(n\) 次的，远大于 \(k\)，所以可以考虑设 \(f(x) = x^n\)，然后将 \(f(x)\) 表示成 \(p_A(x) q(x) + r(x)\) 的形式。惊人的发现将 \(A\) 代入之后，\(f(A) = \mathbf{0} q(A) + r(A) = r(A)\)！也就是说，我们只需要求出来 \(r(x)\) 的多项式然后将 \(A\) 代入即可。

求特征多项式是 \(\mathcal{O}(k^3)\) 或 \(\mathcal{O}(k^4)\) 的，求出 \(r(x)\) 是 \(\mathcal{O}(k^2 \log n)\) 的（直接暴力多项式乘法/多项式取模），将 \(A\) 代回 \(r(x)\) 是 \(\mathcal{O}(k^4)\) 的，总复杂度是 \(\mathcal{O}(k^4 + k^2 \log n)\)。

核心代码：

struct Poly{
	vector<int>v;
	Poly(){v.clear();}
	Poly(vector<int>V){v=V;}
	friend Poly operator*(Poly x,Poly y);
	friend Poly operator%(Poly x,Poly y);//多项式乘法/取模 
}p,x({1}),r({0,1});//p 即 p_A(x) 
int det(int b[][55],int n);//求 n 阶方阵 b 的行列式 
int ans[55];
void gauss(int b[][55],int n);//高斯消元 
int b[55][55],B[55][55];
void gpoly(int a[][55],int n){//暴力求行列式，暴力高斯消元求出 p_A(x) 的各项系数 
	fo(i,0,n){
		fo(j,1,n)fo(k,1,n)b[j][k]=fminus((j==k?i:0),a[j][k]);
		B[i+1][n+2]=det(b,n);
		fo(j,1,n+1)B[i+1][j]=qp(i,j-1);
	}
	gauss(B,n+1),p.v.assign(n+1,0);
	fo(i,0,n)p.v[i]=ans[i+1];
}
string n;
int k,a[55][55],tmp[2][55][55],Ans[55][55];
void solve(){
	cin>>n>>k;
	fo(i,1,k)fo(j,1,k)cin>>a[i][j];
	gpoly(a,k);
	fd(i,n.size()-1,0){//求 r(x) 
		if(n[i]-48)x=x*r%p;
		r=r*r%p;
	}
	fo(i,1,k)tmp[0][i][i]=1;
	fo(i,0,x.v.size()-1){//将 A 代回 r(x) 
		fo(j,1,k)fo(K,1,k)cplus(Ans[j][K],tmp[i&1][j][K]*x.v[i]%mod);
		mem(tmp[!(i&1)],0);
		fo(I,1,k)fo(j,1,k)fo(K,1,k)cplus(tmp[!(i&1)][I][j],tmp[i&1][I][K]*a[K][j]%mod);
	}
	fo(i,1,k)fo(j,1,k)cout<<Ans[i][j]<<" \n"[j==k];
}

关于 Cayley–Hamilton 定理的证明

下文中的 \(n\) 均指方阵 \(A\) 的维度。

设 \(p_A(\lambda) = a_0 \lambda^n + a_1 \lambda^{n - 1} + \dots + a_{n - 1} \lambda + a_n\)，再设 \(B = \lambda E - A\)，则 \(p_A = |B|\)。再设 \(B\) 的伴随矩阵 \(B^* = b_0 \lambda^{n - 1} + b_1 \lambda^{n - 2} + \dots + b_{n - 2} \lambda + b_{n - 1}\)。

由 \(B B^* = |B| E\)，有 \((\lambda E - A)(b_0 \lambda^{n - 1} + b_1 \lambda^{n - 2} + \dots + b_{n - 2} \lambda + b_{n - 1}) = (a_0 \lambda^n + a_1 \lambda^{n - 1} + \dots + a_{n - 1} \lambda + a_n)E\)。

将左式展开：\(b_0 \lambda^n + (b_1 - A b_0) \lambda^{n - 1} + \dots + (b_{n - 1} - A b_{n - 2}) \lambda - A b_{n - 1}\)。则有：

\[\begin{aligned} b_0 = & a_0 E \\ b_1 - A b_0 = & a_1 E \\ \vdots & \\ b_{n - 1} - A b_{n - 2} = & a_{n - 1} E \\ A b_{n - 1} = & a_n E \end{aligned} \]

将第一个式子等号两边乘于 \(A^n\)，第二个式子等号两边乘于 \(A^{n - 1}\)……然后再把左右两边全加起来，可以发现左边加起来全部抵消了，结果是 \(\mathbf{0}\)，而右边是 \(a_0 A^n + a_1 A^{n - 1} + \dots + a_{n - 1} A + a_n = p_A(A)\)，则 \(p_A(A) = \mathbf{0}\)，证毕。