poj3233Matrix Power Series(矩阵乘法)
Matrix Power Series
Time Limit: 3000MS | Memory Limit: 131072K | |
Total Submissions: 23187 | Accepted: 9662 |
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
/* 矩阵乘法经典+二分 Sk=A+A2+A3+...+Ak =(1+Ak/2)*(A+A2+A3+...+Ak/2)+{Ak} =(1+Ak/2)*(Sk/2)+{Ak} 当k为偶数时不要大括号里面的数 */ #include<iostream> #include<cstdio> #include<cstring> using namespace std; int n,m,k; struct matrix { int a[30][30]; void init() { memset(a,0,sizeof a); for(int i=0;i<30;i++) a[i][i]=1; } }; void print(matrix s) { for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { if (j) printf(" "); printf("%d",s.a[i][j]%m); } printf("\n"); } } matrix m_add(matrix a,matrix b)//加法 { matrix c; for(int i=0;i<n;i++) for(int j=0;j<n;j++) c.a[i][j]=((a.a[i][j]+b.a[i][j])%m); return c; } matrix m_mul(matrix a,matrix b)//乘法 { matrix c; for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { c.a[i][j]=0; for(int k=0;k<n;k++) c.a[i][j]+=((a.a[i][k]*b.a[k][j])%m); c.a[i][j]%=m; } } return c; } matrix mul(matrix s,int k)//矩阵快速幂 { matrix ans;ans.init(); while(k>=1) { if(k&1) ans=m_mul(ans,s); k>>=1; s=m_mul(s,s); } return ans; } matrix sum(matrix s,int k)//矩阵前k项求和 { if(k==1) return s; matrix tmp;tmp.init(); tmp=m_add(tmp,mul(s,k>>1));//计算1+A^(k/2) tmp=m_mul(tmp,sum(s,k>>1));//计算(1+A^(k/2))*(A+A^2+A^3+...+A^(k/2)) if(k&1) tmp=m_add(tmp,mul(s,k)); return tmp; } int main() { while(cin>>n>>k>>m) { matrix s; for(int i=0;i<n;i++) for(int j=0;j<n;j++) scanf("%d",&s.a[i][j]); s=sum(s,k); print(s); } }
折花枝,恨花枝,准拟花开人共卮,开时人去时。
怕相思,已相思,轮到相思没处辞,眉间露一丝。