POJ 3233 Matrix Power Series (矩阵+二分+二分)
题目地址:http://poj.org/problem?id=3233
题意:给你一个矩阵A,让你求A+A^2+……+A^k模p的矩阵值
题解:我们知道求A^n我们可以用二分-矩阵快速幂来求,而
当k是奇数A+A^2+……+A^k=A^(k/2+1)+(A+A^2+……A^(k/2))*(1+A^(k/2+1))
当k是偶数A+A^2+……+A^k=(A+A^2+……A^(k/2))*(1+A^(k/2))
可以在一次用二分。
AC代码:
#include <iostream> #include <cstdio> #include <cstring> #include <string> #include <cstdlib> #include <cmath> #include <vector> #include <list> #include <deque> #include <queue> #include <iterator> #include <stack> #include <map> #include <set> #include <algorithm> #include <cctype> using namespace std; typedef long long LL; const int N=31; const int mod=1000007; const int INF=0x3f3f3f3f; const double PI=acos(-1.0); int n,k,p; struct M { int m[N][N]; }; void print(M t) { int i,j; for(i=1;i<=n;i++) { for(j=1;j<n;j++) printf("%d ",t.m[i][j]); printf("%d\n",t.m[i][n]); } } M xh_mod(M a) { M t; int i,j; for(i=1;i<=n;i++) for(j=1;j<=n;j++) t.m[i][j]=a.m[i][j]%p; return t; } M xh_mult(M a,M b) { M t; int i,j,k; memset(t.m,0,sizeof(t.m)); for(i=1;i<=n;i++) for(j=1;j<=n;j++) for(k=1;k<=n;k++) t.m[i][j]=(t.m[i][j]+a.m[i][k]*b.m[k][j])%p; return t; } M xh_pow(M a,int b) { M t; memset(t.m,0,sizeof(t.m)); for(int i=1;i<=n;i++) t.m[i][i]=1; while(b) { if(b&1) t=xh_mult(t,a); a=xh_mult(a,a); b/=2; } return t; } M xh_add(M a,M b) { M t; int i,j; for(i=1;i<=n;i++) for(j=1;j<=n;j++) t.m[i][j]=(a.m[i][j]+b.m[i][j])%p; return t; } M love(M a,int k) { M t,x; if(k==1) { t=a; return t; } x=love(a,k/2); if(k&1) { M o=xh_pow(a,k/2+1); return xh_add(xh_add(x,o),xh_mult(x,o)); } else { M o=xh_pow(a,k/2); return xh_add(x,xh_mult(x,o)); } } int main() { int i,j; while(~scanf("%d%d%d",&n,&k,&p)) { M a,t; for(i=1;i<=n;i++) for(j=1;j<=n;j++) scanf("%d",&a.m[i][j]); t=xh_mod(a); a=love(t,k); print(a); } return 0; }