poj 3233 Matrix Power Series 矩阵构造+快速幂

点击打开链接

要求A+A2+...+Ak，而不是单个矩阵的幂

即求A(E+A+..A^k-1) 可以不断乘以K次A

那么我们可以构造一个分块的辅助矩阵 S，其中 A 为原矩阵，E 为单位矩阵，O 为0矩阵

即 S=[A,E]   容易用归纳得到S^k=[A^k,A^k-1+A^k-2+...A+E]
 [0 ,E]
[0,E]
#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
using namespace std;
typedef long long ll;
const int N=80;
int a

,b

,n,k,mod;
void mul(int x

,int y

)
{
int i,j,k;
int c

;
for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
{
c[i][j]=0;
for(int k=1;k<=n;k++)
{
c[i][j]=(c[i][j]+x[i][k]*y[k][j])%mod;
}
}
}
for (i=1;i<=n;i++)
for (j=1;j<=n;j++)
x[i][j]=c[i][j];
}
void quickpow(int n)
{
while(n)
{
if(n&1)
mul(b,a);

mul(a,a);
n>>=1;
}
}
int main()
{
cin>>n>>k>>mod;
memset(a,0,sizeof(a));
for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
{
cin>>a[i][j];
a[i][i+n]=a[i+n][i+n]=b[i][i]=b[i+n][i+n]=1;
//b为单位阵
}
}//构造分块矩阵
n*=2;
quickpow(k+1);

n/=2;
for(int i=1;i<=n;i++)
b[i][i+n]--;

for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
{
printf("%d",b[i][j+n]);
if(j==n)
cout<<endl;
else
cout<<' ';
}
}
return 0;
}