描述
给定n × n矩阵A和正整数k,求和 。
输入
输入只包含一个测试用例。输入的第一行包含三个正整数N(N ≤30),k(k ≤)和m(m <
)。然后跟随n行,每行包含32,768以下的n个非负整数,以行 - 主顺序给出A的元素。
产量
以与给定A相同的方式输出S模数m的元素。
样本输入
2 2 4
0 1
1 1
样本输出
1 2
2 3
解法1:矩阵快速幂+二分
题目大意:给定矩阵A,求A + A^2 + A^3 + … + A^k的结果(两个矩阵相加就是对应位置分别相加)。输出的数据mod m。k<=10^9。
这道题两次二分,相当经典。首先我们知道,A^i可以二分求出。然后我们需要对整个题目的数据规模k进行二分。比如,当k=6时,有:
A + A^2 + A^3 + A^4 + A^5 + A^6 =(A + A^2 + A^3) + A^3*(A + A^2 + A^3)
应用这个式子后,规模k减小了一半。我们二分求出A^3后再递归地计算A + A^2 + A^3,即可得到原问题的答案。
#include <iostream>
#include <cstdio>
#include <string>
#include <cstring>
#include <stdlib.h>
#include <math.h>
#include <ctype.h>
#include <queue>
#include <map>
#include <set>
#include <algorithm>
using namespace std;
int mod;
int n;
struct matrix
{
int ma[40][40];
}init, res;
matrix Mult(matrix x, matrix y)
{
int i, j, k;
matrix tmp;
memset(tmp.ma,0,sizeof(tmp.ma));
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
for(k=0;k<n;k++)
{
tmp.ma[i][j]=(tmp.ma[i][j]+x.ma[i][k]*y.ma[k][j])%mod;
}
}
}
return tmp;
}
matrix Pow(matrix x, int k)
{
matrix tmp;
int i, j;
memset(tmp.ma,0,sizeof(tmp.ma));
for(i=0;i<n;i++) tmp.ma[i][i]=1;
while(k)
{
if(k&1) tmp=Mult(tmp,x);
x=Mult(x,x);
k>>=1;
}
return tmp;
}
matrix add(matrix x, matrix y)
{
int i, j;
matrix tmp;
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
tmp.ma[i][j]=(x.ma[i][j]+y.ma[i][j])%mod;
}
}
return tmp;
}
matrix sum(matrix x, int k)
{
matrix tmp, y;
if(k==1) return x;
tmp=sum(x,k/2);
if(k&1)
{
y=Pow(x,k/2+1);
tmp=add(Mult(y,tmp),tmp);
return add(tmp,y);
}
else
{
y=Pow(x,k/2);
return add(Mult(y,tmp),tmp);
}
}
int main()
{
int k, m, x, i, j;
scanf("%d%d%d",&n,&k,&mod);
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
scanf("%d",&x);
init.ma[i][j]=x%mod;
}
}
res=sum(init, k);
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
printf("%d",res.ma[i][j]);
if(j!=n-1) printf(" ");
}
puts("");
}
return 0;
}
解法二:构造矩阵+矩阵快速幂
设S[k]=A + A2 + A3 + … + Ak.
那么S[n]=A*S[n-1]+A;
=>所求 S=(A,1)*(A,0 )^(k-1)
A,1
注意:矩阵套矩阵,里面的1代表单位矩阵,即对角线上元素为1,其余为0。
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <vector>
#include <queue>
#include <algorithm>
using namespace std;
typedef long long ll;
typedef struct {
ll m[120][120];
}Matrix;
ll n,k,mod;
Matrix Mul(Matrix a, Matrix b)
{
Matrix c;
memset(c.m, 0, sizeof(c.m));
for (int i = 0; i < 2*n; i++)
{
for (int j = 0; j < 2*n; j++)
{
for (int k = 0; k < 2*n; k++)
{
c.m[i][j] = (c.m[i][j] + (a.m[i][k] * b.m[k][j]) % mod ) % mod;
}
}
}
return c;
}
//矩阵乘法
Matrix solve(Matrix a, Matrix b)
{
Matrix c;
memset(c.m, 0, sizeof(c.m));
for (int i = 0; i < n; i++)
{
for (int j = 0; j < 2*n; j++)
{
for (int k = 0; k < 2*n; k++)
{
c.m[i][j] = (c.m[i][j] + (a.m[i][k] * b.m[k][j]) % mod ) % mod;
}
}
}
return c;
}
Matrix fastm(Matrix a, ll num)
{
Matrix res;
memset(res.m, 0, sizeof(res.m));
//初始化为单位矩阵
for(int i=0;i<2*n;i++)
res.m[i][i]=1;
while (num)
{
if (num & 1)
res = Mul(res, a);
num >>= 1;
a = Mul(a, a);
}
return res;
}
int main()
{
Matrix a;
memset(a.m,0,sizeof(a.m));
scanf("%lld%lld%lld",&n,&k,&mod);
for(int i=0;i<n;i++)
{
for(int j=0;j<n;j++)
{
scanf("%lld",&a.m[i][j]);
a.m[i][j]=a.m[i][j]%mod;
}
}
Matrix b;
memset(b.m,0,sizeof(b.m));
for(int i=0;i<n;i++)
for(int j=0;j<n;j++)
b.m[i][j]=a.m[i][j];
for(int i=0;i<n;i++)
for(int j=n;j<2*n;j++)
b.m[i][j]=0;
for(int i=n;i<2*n;i++)
for(int j=0;j<n;j++)
b.m[i][j]=a.m[i-n][j];
//注意:矩阵等于1表示该矩阵是单位矩阵,即主对角线为1,其余为0
for(int i=n;i<2*n;i++)
b.m[i][i]=1;
for(int i=0;i<n;i++)
a.m[i][i+n]=1;
b=fastm(b,k-1);
Matrix ans=solve(a,b);
for(int i=0;i<n;i++)
{
printf("%lld",ans.m[i][0]);
for(int j=1;j<n;j++)
{
printf(" %lld",ans.m[i][j]);
}
printf("\n");
}
return 0;
}