POJ3233 Matrix Power Series（矩阵）

最新推荐文章于 2020-10-17 09:49:19 发布

Xuh7

最新推荐文章于 2020-10-17 09:49:19 发布

阅读量509

点赞数

分类专栏： # 动态规划算法文章标签： poj 矩阵快速幂

本文链接：https://blog.csdn.net/xuh723/article/details/40952601

版权

算法同时被 2 个专栏收录

31 篇文章 0 订阅

订阅专栏

动态规划

3 篇文章 0 订阅

订阅专栏

题目链接：http://poj.org/problem?id=3233

Matrix Power Series

Time Limit: 3000MS		Memory Limit: 131072K
Total Submissions: 15856		Accepted: 6782

Description

Given a n × n matrix A and a positive integer k, find the sum S = A + A² + A³ + … + A^k.

Input

The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 10⁹) and m (m < 10⁴). Then follow n lines each containing n nonnegative integers below 32,768, giving A’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和m，求S=A+A^2+A^3+...+A^k，输出S模m的值。

矩阵的快速幂时间复杂度为O(n3logk)，但是由于还要累加，因此复杂度变大。

事实上，只需要把原来的矩阵补成一个分块矩阵，即可快速求出累加的和。

分块为 A 0

I I

这个剧真的k次幂乘以分块矩阵(I 0)，即为分块矩阵（A^k,Sk），其中Sk=I+A+A^2+...+A^k

减去I即为S。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
using namespace std;

int a[35][35];
int b[70][70],ans[70][70];
int n,k,m;

void Qmod(int x)
{
	memset(ans,0,sizeof(ans));
	for (int i=0;i<2*n;i++) ans[i][i]=1;
	while (x)
	{
		if (x&1)
		{
			int tem[70][70];
			for (int i=0;i<2*n;i++)
			{
				for (int j=0;j<2*n;j++)
				{
					tem[i][j]=0;
					for (int k=0;k<2*n;k++)
					{
						tem[i][j]+=ans[i][k]*b[k][j];
						tem[i][j]%=m;
					}
				}
			}
			for (int i=0;i<2*n;i++)
				for (int j=0;j<2*n;j++)
					ans[i][j]=tem[i][j];
		}
		int tem[70][70];
		for (int i=0;i<2*n;i++)
		{
			for (int j=0;j<2*n;j++)
			{
				tem[i][j]=0;
				for (int k=0;k<2*n;k++)
				{
					tem[i][j]+=b[i][k]*b[k][j];
					tem[i][j]%=m;
				}
			}
		}
		for (int i=0;i<2*n;i++)
			for (int j=0;j<2*n;j++)
				b[i][j]=tem[i][j];
		x>>=1;
	}
}

int main()
{
	while (scanf("%d%d%d",&n,&k,&m)!=EOF)
	{
		memset(b,0,sizeof(b));
		for (int i=0;i<n;i++)
		{
			for (int j=0;j<n;j++)
			{
				scanf("%d",&a[i][j]);
				b[i][j]=a[i][j];
			}
		}
		for (int i=0;i<n;i++)
		{
			b[i+n][i]=1;b[i+n][i+n]=1;//制造分块矩阵
		}
		Qmod(k+1);
		for (int i=0;i<n;i++)
		{
			for (int j=0;j<n-1;j++)
			{
				if (i==j) ans[i+n][j]=(ans[i+n][j]-1+m)%m;//减去单位矩阵I
				cout<<ans[i+n][j]<<" ";
			}
			if (i==n-1) ans[i+n][n-1]=(ans[i+n][n-1]-1+m)%m;
			cout<<ans[i+n][n-1]<<endl;
		}
	}
	return 0;
}