題目鏈接：http://poj.org/problem?id=3233

題目意思：給一個矩陣n*n的矩陣A和一個k，求一個式子 S = A + A² + A³ + … + A^k。

這個需要用到等比數列和的二分加速。

當n為奇數的時候，Sn=Sn-1+A^k;

當n為偶數的時候，Sn=（S[n/2]+E）*A^(k/2)

自己xjb推一下就知道等比數列和的二分加速是咋回事了。我舉個例子，我們假設求等比數列2,4,8,16,32,64的和s=（8+1）*（2+4+8），而2+4+8=（2+4）+8，而（2+4）=(2+1)*2,其他的可以一次類推。

這個過程我們可以直接用一個遞歸過程算出來，其余我們套用矩陣快速冪的模板就好了。

代碼：

 1 //Author: xiaowuga
 2 #include<iostream>
 3 #include<cstring>
 4 #define maxx INT_MAX
 5 #define minn INT_MIN
 6 #define inf 0x3f3f3f3f
 7 #define size 35 
 8 using namespace std;
 9 typedef long long ll;
10 int n,k,mod;
11 struct Mat{
12     int mat[size][size];
13     void clear(){
 
14         memset(mat,0,sizeof(mat));
15     }
16 
17     Mat operator *(const Mat &e) const{
18         Mat tmp;
19         tmp.clear();
20         for(int k=0;k<n;k++)
21             for(int i=0;i<n;i++){
22                 if(mat[i][k]==0) continue;
23                 for(int j=0;j<n;j++){
 
24                    if(e.mat[k][j]==0) continue;
25                     tmp.mat[i][j]+=mat[i][k]*e.mat[k][j]%mod;
26                     tmp.mat[i][j]%=mod;
27                 }
28             }
29         return tmp;
30     }
31     Mat operator +(const Mat &e) const{
32         Mat tmp;
33         tmp.clear();
34             for(int i=0;i<n;i++){
35                 for(int j=0;j<n;j++){
36                     tmp.mat[i][j]=(mat[i][j]%mod+e.mat[i][j]%mod)%mod;
37                 }
38             }
39         return tmp;
40     }
41 };
42 Mat m,E;
43 Mat pow(Mat ma,ll num){
44     Mat ans;
45     ans.clear();
46     for(int i=0;i<n;i++) ans.mat[i][i]=1;
47     while(num){
48         if(num&1) ans=ans*ma;
49         num/=2;
50         ma=ma*ma;
51     }
52     return ans;
53 }
54 Mat fun(int x){
55     if(x==1) return m;
56     if(x&1) return fun(x-1)+pow(m,x);
57     else return (pow(m,x/2)+E)*fun(x/2); 
58 }
59 int main() {
60     ios::sync_with_stdio(false);cin.tie(0);
61     cin>>n>>k>>mod;
62     E.clear();
63     for(int i=0;i<n;i++) E.mat[i][i]=1;
64     for(int i=0;i<n;i++)
65         for(int j=0;j<n;j++) cin>>m.mat[i][j];
66     Mat ans=fun(k);
67     for(int i=0;i<n;i++){
68         for(int j=0;j<n;j++)
69             cout<<ans.mat[i][j]<<" ";
70         cout<<endl;
71     }
72     return 0;
73 }

poj3233—Matrix Power Series