矩阵儿快速幂 - POJ 3233 矩阵力量系列

不要管上面的标题的bug
那是幂的意思,不是力量。。。


POJ 3233 Matrix Power Series

描述

Given a n × n matrix A and a positive integer k, find the sum $ S = A + A^2 + A^3 + … + A^k $.
给你个n×n大小的矩阵A和一个正整数k,求矩阵S = A + A^2 + A^3 + … + A^k。

输入

The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 109) and m (m < 104).
输入只包含一个测试点。 第一行输入包含三个正整数n(n≤30),k(k≤10^9)和m(m <10^4)。
Then follow n lines each containing n nonnegative integers below 32,768, giving A’s elements in row-major order.
接下来n行,每行有n个低于32,768的非负整数,描述A的组成。

输出

Output the elements of S modulo m in the same way as A is given.
把矩阵S的值按输入A的格式模上m输出。

样例

.in
2 2 4
0 1
1 1
.out
1 2
2 3

暴力求解的话要交给天河2号评测
吾等平民还是写点高效算法吧
先把那个式子拆成下面这样子:
$ \sum^k_{i=1} {A^i} = (\sum^{k'}_{i=1} {A^i}) \times {(1 + A^{k'+1})} + 2|k?A^k:0 $
然后DFS分一分合一下就好了。

代码蒯上

#include<iostream>
#include<iomanip>
#include<cstdio>
#include<cstdlib>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;
inline int gotcha()
{
	register int a=0;bool b=1;register char c=getchar();
	while(c>'9' || c<'0'){if(c=='-')b=0;c=getchar();}
	while(c>='0' && c<='9')a=a*10+c-48,c=getchar();
	return b?a:-a;
}
int n,tim,mo;
struct martix
{
	int a[32][32];
	martix(){memset(a,0,sizeof(a));}
	const martix operator + (const martix &b)const
	{
		martix c;register int i,j;
		for(i=1;i<=n;i++)for(j=1;j<=n;j++)
			c.a[i][j]=a[i][j]+b.a[i][j],c.a[i][j]%=mo;
		return c;
	}
	const martix operator * (const martix &b)const
	{
		martix c;register int i,j,k;
		for(k=1;k<=n;k++)for(i=1;i<=n;i++)for(j=1;j<=n;j++)
			c.a[i][j]+=a[i][k]*b.a[k][j],c.a[i][j]%=mo;
		return c;
	}
}ori;
martix powa(martix in,int tim)
{
	if(tim==1)return in;
	martix tmp;
	for(int i=0;i<=n;i++)tmp.a[i][i]=1;
	if(tim==0)return tmp;
	while(tim){if(tim&1)tmp=in*tmp;in=in*in;tim>>=1;}
	return tmp;
}
martix finder(martix now,int k)
{
	if(k==1)return now;
	if(k&1)return finder(now,k-1)+powa(now,k);
	return (powa(now,0)+powa(now,k>>1))*finder(now,k>>1);
}
int main()
{
	register int i,j;
	n=gotcha(),tim=gotcha(),mo=gotcha();
	for(i=1;i<=n;i++)for(j=1;j<=n;j++)ori.a[j][i]=gotcha()%mo;
	ori=finder(ori,tim);
	for(i=1;i<=n;i++){for(j=1;j<=n;j++)printf("%d ",ori.a[j][i]);printf("\n");}
	return 0;
}
posted @ 2017-10-10 21:27  iot;  阅读(169)  评论(2编辑  收藏  举报
知识共享许可协议
年轻人,你需要更多的知识