《ybtoj高效进阶》第六部分第一章例题4 矩阵求和

题目大意

给方形(长为n)矩阵A，求($A+A^2+A^3+\dots +A^k)\mathrm{m}\mathrm{o}\mathrm{d}m$

思路

设$B(2n\ast 2n)$左上角为A,右上，右下角为单位矩阵，剩下为0，显然我们的结果是$B^{k+1}$的右上角减一个单位矩阵。
code：

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
using namespace std;
long long n,myd=1e9+7,m; 
const int lim=101;
struct Matrix{
	long long m[lim+1][lim+1];
	bool wj[lim+1][lim+1];
} k;
Matrix cf(Matrix &x,Matrix &y,long long myd)
{
	Matrix z;
	memset(z.m,0,sizeof(z.m));
	for (int i=0;i<=lim;i++) for (int j=0;j<=lim;j++) if (x.m[i][j])
	{
		for (int k=0;k<=lim;k++)
		{
			z.wj[i][k]|=((z.m[i][k]+(x.m[i][j]*y.m[j][k]))>=myd)|(x.wj[i][j])|(y.wj[i][j]);
			z.m[i][k]=(z.m[i][k]+(x.m[i][j]*y.m[j][k])%myd)%myd;
		}
	}
	return z;
}
Matrix jzksm(Matrix x,long long y,long long myd)
{
	Matrix ans;
	memset(ans.m,0,sizeof(ans.m));
	for (int i=0;i<=lim;i++) ans.m[i][i]=1;
	while (y)
	{
		if (y&1) ans=cf(ans,x,myd);
		x=cf(x,x,myd);
		y>>=1;
	}
	return ans;
}
Matrix a;
int main()
{
	cin>>n>>m>>myd;
	for (int i=0;i<n;i++) for (int j=0;j<n;j++) cin>>a.m[i][j];
	for (int i=n;i<2*n;i++) a.m[i][i]=1;
	for (int i=0;i<n;i++) a.m[i][i+n]=1;
	k=jzksm(a,m+1,myd);
	for (int i=0;i<n;i++)
	{
		for (int j=n;j<n*2;j++)
		{
			if (j-i==n) k.m[i][j]=k.m[i][j]-1;
			cout<<k.m[i][j]<<' ';
		}
		cout<<endl;
	}
	return 0;
}