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
",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;
}