POJ 3233 Matrix Power Series(マトリクス+二分+二分)

タイトルアドレス:http://poj.org/problem?id=3233
标题:A+A^2+……+A^kモードpのマトリクス値を求めるマトリクスAをあげる
A^nを求めるには二分行列の高速べき乗で求めることができることを知っています.
kが奇数A+A^2+…+A^k=A^(k/2+1)+(A+A^2+…(k/2))*(1+A^(k/2+1))
kが偶数A+A^2+…+A^k=(A+A^2+…...A^(k/2)*(1+A^(k/2))
一度に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;

}