[caioj 1481及poj 3233，利用矩阵乘法解决的经典题目三]矩阵幂级数

给定一个n行n列的矩阵A，求A+A^2+A^3+…+A^k的结果，并且输出的每个数都mod m。

这道题非常经典，因为它需要两次快速幂，为什么？首先A^i我们当然可以求出来，可以用快速幂，但如果计算k次就会超时。所以我们要再进行一次，举个例子，A+A^2+A^3+…+A^6，可以转化成A+A^2+A^3+A^3(A+A^2+A^3),这样就可以先二分算出(A+A^2+A^3)，然后再乘A^3，便可以避免计算(A^4+A^5+A^6)，速度就可以大大提升。如果k是奇数，只需再加上A^k即可。

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
struct node
{
    int v[35][35];
    node()
    {
        memset(v,0,sizeof(v));
    }
};
int n,m;node per,a;
node jiafa(node a,node b)//直接各个数相加即可
{
    node c;
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++)
        {
            c.v[i][j]=(c.v[i][j]+a.v[i][j]+b.v[i][j])%m;
        }
    }
    return c;
}
node chengfa(node a,node b)
{
    node c;
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++)
        {
            for(int k=1;k<=n;k++)
            {
                c.v[i][j]=(c.v[i][j]+a.v[i][k]*b.v[k][j])%m;
            }
        }
    }
    return c;
}
node power(int n)//第一次快速幂
{
    node ans=per;
    node A=a;
    while(n>0)
    {
        if(n%2==1)ans=chengfa(ans,A);
        A=chengfa(A,A);
        n=n/2;
    }
    return ans;
}
node getsum(int k)//第二次快速幂
{
    if(k==1)return a;
    node b,t;
    if(k%2==0)
    {
        t=getsum(k/2);
        t=jiafa(t,chengfa(t,power(k/2)));
    }
    else
    {
        t=getsum(k/2);
        t=jiafa(t,chengfa(t,power(k/2)));
        t=jiafa(t,power(k));
    }
    return t;
}
int main()
{
    int i,j,k,x;
    scanf("%d%d%d",&n,&k,&m);   
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++)
        {
            scanf("%d",&a.v[i][j]);
            a.v[i][j]%=m;
            if(i==j)per.v[i][j]=1;//构建单位矩阵
            else per.v[i][j]=0;
        }
    }
    node ans=getsum(k);
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<n;j++)printf("%d ",ans.v[i][j]);
        printf("%d\n",ans.v[i][n]);
    }
    return 0;
}