POJ-3233 Matrix Power Series（矩阵快速幂）

Matrix Power Series

Time Limit: 3000MS		Memory Limit: 131072K
Total Submissions: 23032		Accepted: 9606

Description

Given a n × n matrix A and a positive integer k, find the sum S = A + A² + A³ + … + A^k.

Input

The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 10⁹) and m (m < 10⁴). 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

题解：简单矩阵快速幂

把A看成一个数，可以构建这样一个矩阵

|A 0|    *    |  A  |

|E E|         |sum|

其中sum为所求答案，E为单位矩阵

#include<cstdio>
#include<algorithm>
#include<string.h>
using namespace std;
int n,mod;
struct Mat{
    int x[31][31],sz;
    Mat(){
        memset(x,0,sizeof(x));
        sz=n;
    }
    void init(){
        for(int i=0;i<sz;i++) x[i][i]=1;
    }
    Mat operator*(const Mat& m)const{
        Mat ret;
        for(int i=0;i<sz;i++)
            for(int j=0;j<sz;j++)
                for(int k=0;k<sz;k++)
                    ret.x[i][j]=(ret.x[i][j]+x[i][k]*m.x[k][j])%mod;
        return ret;
    }
    Mat operator+(const Mat& m)const{
        Mat ret;
        for(int i=0;i<sz;i++)
            for(int j=0;j<sz;j++)
                ret.x[i][j]=(ret.x[i][j]+x[i][j]+m.x[i][j])%mod;
        return ret;
    }
    void print(){
        for(int i=0;i<sz;i++)
            for(int j=0;j<sz;j++)
                printf("%d%c",x[i][j],j==sz-1?'\n':' ');
    }
};
struct MMat{
    Mat a[2][2];
    MMat(){
        for(int i=0;i<2;i++)
            for(int j=0;j<2;j++)
                memset(a[i][j].x,0,sizeof(a[i][j].x));
    }
    void init(){
        for(int i=0;i<2;i++) a[i][i].init();
    }
    MMat operator*(const MMat& m)const{
        MMat ret;
        for(int i=0;i<2;i++)
            for(int j=0;j<2;j++)
                for(int k=0;k<2;k++)
                    ret.a[i][j]=ret.a[i][j]+a[i][k]*m.a[k][j];
        return ret;
    }
};
MMat pow(MMat A,int k){
    MMat ret;
    ret.init();
    while(k){
        if(k&1) ret=ret*A;
        A=A*A;
        k>>=1;
    }
    return ret;
}
Mat solve(int k){
    Mat A;
    for(int i=0;i<n;i++)
        for(int j=0;j<n;j++)
            scanf("%d",&A.x[i][j]);
    MMat MA,MB;
    MA.a[0][0]=A;
    MA.a[1][0].init();MA.a[1][1].init();
    MB.a[0][0]=A;
    MA=pow(MA,k);
    return (MA*MB).a[1][0];
}
int main(){
    int k;
    //freopen("in.txt","r",stdin);
    while(~scanf("%d%d%d",&n,&k,&mod)){
        Mat A=solve(k);
        A.print();
    }
    return 0;
}