Matrix Power Series
Time Limit:3000MS Memory Limit:131072K
Description
Given a n × n matrix A and a positive integer k, find the sum S = A + A2 + A3 + … + Ak.
Input
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). 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
Source
POJ Monthly--2007.06.03, Huang, Jinsong
题意：给出矩阵A，求S=A+A2+A3+...+Ak的值。
题解：既然是要求矩阵多个幂之和，那么不用说肯定得用矩阵快速幂了。
但是针对此题，如果暴力计算每一个矩阵的幂，很显然是要TLE的，因为k的规模高达109。
于是，我们便想到了用二分的方法优化。
对于S，我们可以进行如下拆分：
S(k)=A+A2+A3+...+Ak=(1+Ak/2)×(A+A2+A3+...+Ak/2)（k为偶数）
S(k)=A+A2+A3+...+Ak=(1+Ak/2)×(A+A2+A3+...+Ak/2)+Ak（k为奇数，其中除号表示相除后向下取整）
我们发现拆项中含有A+A2+A3+...+Ak/2，即S(k/2)，于是只要不断二分即可求解。
1 #include<cstdio> 2 #include<cstring> 3 using namespace std; 4 const int N=30; 5 int n,k,mod; 6 struct node{ 7 int a[N][N]; 8 void clr(){memset(a,0,sizeof(a));} 9 void init(){for (int i=0;i<n;++i) a[i][i]=1;} 10 }base,ans; 11 node operator + (node a,node b) 12 { 13 node c; c.clr(); 14 for (int i=0;i<n;++i) 15 for (int j=0;j<n;++j) 16 c.a[i][j]=(a.a[i][j]+b.a[i][j])%mod; 17 return c; 18 } 19 node operator * (node a,node b) 20 { 21 node c; c.clr(); 22 for (int i=0;i<n;++i) 23 for (int j=0;j<n;++j) 24 for (int k=0;k<n;++k) 25 c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%mod; 26 return c; 27 } 28 node Pow(node now,int x) 29 { 30 if (x==1) return now; 31 node tmp; --x; 32 for (tmp=now;x;x>>=1,now=now*now) 33 if (x&1) tmp=tmp*now; 34 return tmp; 35 } 36 node S(node now,int x) 37 { 38 if (x==1) return now; 39 node tmp; tmp.clr(); tmp.init(); 40 tmp=tmp+Pow(now,x/2); 41 tmp=tmp*S(now,x/2); 42 if (x%2) tmp=tmp+Pow(now,x); 43 return tmp; 44 } 45 int main() 46 { 47 scanf("%d%d%d",&n,&k,&mod); 48 base.clr(); 49 for (int i=0;i<n;++i) 50 for (int j=0;j<n;++j) 51 scanf("%d",&base.a[i][j]); 52 ans=S(base,k); 53 for (int i=0;i<n;++i) 54 { 55 for (int j=0;j<n;++j) 56 printf("%d ",ans.a[i][j]); 57 printf(" "); 58 } 59 return 0; 60 }