POJ 3233 Matrix Power Series

Description

Given a n × n matrix A and a positive integer k, find the sumS = A + A2 + A3 + … + Ak.

Input

The input contains exactly one test case. The first line of input contains three positive integersn (n ≤ 30), k (k ≤ 109) and m (m < 104). Then follown lines each containing n nonnegative integers below 32,768, givingA’s elements in row-major order.

Output

Output the elements of S modulo m in the same way as A is given.

Sample Input

2 2 4
0 1
1 1

Sample Output

1 2
2 3


题意:给定n*n的矩阵A,整数k,mod值为m。求S=A+A^2+...+A^k。输出S矩阵,并且S的每一个元素对m取余数。


思路:

递推式A+...+A^(2k)=(A+...+A^k)+(A+...+A^k)*A^k

原理:A^2k=(A^k)*(A^k)

如果k为偶数,那么(A+A^2+....A^K) = (A+...+A^K/2)+A^K/2*(A+...+A^K/2)

如果k为奇数,那么(A+A^2+....A^K) = (A+...+A^K/2)+A^K/2*(A+...+A^K/2)+A^k

例子:

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+A^2+A^3=A+A^2+A^2*(A)                     (k为数)

代码如下:

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;
int n,k;
int mod; //n*n阶的矩阵;k是sum=A^1+A^2+.....+A^k;m是被模的元素  
struct Matrix{
	int matrix[50][50];
	Matrix(int a=0){
		memset(matrix,0,sizeof(matrix)); //矩阵清零 
		if(a==1) //单位矩阵
			for(int i=0;i<50;i++)
				matrix[i][i]=1;
	}
}m;
Matrix operator * (Matrix a,Matrix b){ //矩阵乘法
	Matrix res;
	memset(res.matrix,0,sizeof(res.matrix));
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			for(int k=1;k<=n;k++)
				res.matrix[i][j]=(res.matrix[i][j]+(a.matrix[i][k]*b.matrix[k][j])%mod)%mod;
	return res;
}
Matrix operator + (Matrix a,Matrix b){ //矩阵加法
	Matrix res;
	for(int i=1;i<=n;i++)
		for(int j=1;j<=n;j++)
			res.matrix[i][j]=(a.matrix[i][j]+b.matrix[i][j])%mod; 
	return res;
}
Matrix operator ^ (Matrix a,int k){ //矩阵快速幂
	bool flag=false;
	Matrix ans=1;
	while(k){
		if(k&1){
			if(flag)ans=ans*a;
			else ans=a;
			flag=true;
		}
		a=a*a;
		k>>=1;
	}
	return ans;
}
Matrix sum(int k){ //求sum( S(k) )= A + A2 + A3 + …+ Ak
	if(k==1)return m; //递归回去 
	else{
		Matrix tmp=sum(k>>1); //sum( S(k/2) )
		if(k&1){ //k为奇数 
			Matrix tmp2=m^((k>>1)+1); 
			return tmp+tmp2+tmp*tmp2;
		}
		else{
			Matrix tmp2=m^(k>>1);
			return tmp+tmp*tmp2; //sum( S(k/2) )*(tmp2+1)
		}
	}
}
int main(){
	int i,j;
	Matrix ans;
	scanf("%d%d%d",&n,&k,&mod);
	for(i=1;i<=n;i++)
		for(j=1;j<=n;j++){
			scanf("%d",&m.matrix[i][j]);
			m.matrix[i][j]%=mod;
		}
	ans=sum(k);
	for(i=1;i<=n;i++,puts(""))
		for(j=1;j<=n;j++)
			printf("%d ",ans.matrix[i][j]);
	return 0;
}




  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 5
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值