题意分析
给定一个 \(n \times n\) 的矩阵 \(A\) 和正整数 \(k\)，求 \(S=A^1+A^2+\cdots+A^k\)。
解题思路
求 \(A^n\) 要用到矩阵快速幂。但是 \(k \le 10^9\)，求 \(k\) 个幂会超时，所以需要用到分治的策略。
我们知道，矩阵乘法满足乘法分配律，即 \(A \times C + B \times C = (A + B) \times C\)。那么我们就可以二分 \(k\)，\(A^{mid+1}\) 加到 \(A^r\) 就等于 \(A^l\) 加到 \(A^{mid}\) 乘 \(A^{\frac{r-l+1}{2}}\)。当 \((r-l+1)\) 是偶数即 \(l\) 到 \(r\) 之间的项数为偶数时，求 \(A^l\) 加到 \(A^{mid}\) 乘 \(A^{\frac{r-l+1}{2}}\) 再加 \(A^l\) 加到 \(A^{mid}\)，而求 \(A^l\) 加到 \(A^{mid}\) 时可以继续二分，直到 \(l=r\) 为止。当 \((r-l+1)\) 是奇数时，我们把多出来的一项单独求（这里选择第 \(l\) 项），剩下的用二分求。如此一来时间复杂度就会压缩到 \(O(\log k)\)。
写一个简单的递归即可。
代码实现
#include<iostream>
#include<cstring>
#define int long long
using namespace std;
struct Matrix{int mx[110][110];
}a,s,c;
int n,m,k;
Matrix operator *(const Matrix &a,const Matrix &b){//运算符重载为矩阵乘法 memset(c.mx,0,sizeof(c.mx));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){for(int d=1;d<=n;d++){c.mx[i][j]=(c.mx[i][j]+(a.mx[i][d]*b.mx[d][j])%m)%m;}}}return c;
}
Matrix operator +(const Matrix &a,const Matrix &b){//运算符重载为矩阵加法 memset(c.mx,0,sizeof(c.mx));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){c.mx[i][j]=(a.mx[i][j]+b.mx[i][j])%m;}}return c;
}
Matrix ksm(Matrix a,int k){//快速幂 Matrix d,e=a;memset(d.mx,0,sizeof(d.mx));for(int i=1;i<=n;i++) d.mx[i][i]=1;while(k){if(k&1) d=d*e;e=e*e;k>>=1;}return d;
}
Matrix sum(int l,int r){//二分 if(l==r){return ksm(a,l);}if((r-l+1)%2==0){int mid=(l+r)/2;Matrix x=ksm(a,mid+1-l),y=sum(l,mid);return y+y*x;}else{return sum(l,l)+sum(l+1,r);}
}
signed main(){cin>>n>>k>>m;for(int i=1;i<=n;i++) s.mx[i][i]=1;for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){cin>>a.mx[i][j];a.mx[i][j]%=m;}}s=sum(1,k);for(int i=1;i<=n;i++){for(int j=1;j<=n;j++) cout<<s.mx[i][j]<<" ";cout<<endl;}return 0;
}
注意事项
记得随时取余 \(m\)。
