poj 3233 等比矩阵的性质来计算举证的快速幂

poj 3233 Matrix Power Series

  这道题大概的题意给定一个n*n的矩阵，和k，m；

要求计算s=A+A^2+A^3+A^4+……+A^k.

其中n<=30,k<=10^9,m<10^4;

输出S%m；

我的理解和分析：

这道题有两种思路：

第一种思路就是等比矩阵的性质，

|A    1|        |A^n1+A+A^2+A^3+……+A^n-1 |  

|0    1|的n次方就是|0 1       |

不信你可以举个栗子，比如说A=2，n=2 ，那么

|2  1|     |4 1+2|

|0 1|的平方就是|0    1|

假设前边的矩阵就是A， 后边的就是B，那么要求的就是

（1+A+A^2+A^3+……+A^n-1）-1+A^n，对吧，就这样。

那么怎么把A矩阵变成B矩阵呢？

主要就是下边，请你们注意观察：

int n;//n就是这个A矩阵的范围

memet(a.m,0,sizeof a.m);//a.m就是这个A矩阵里边的数组的值

for(int i=1;i<=n;i++)
            {
                a.m[i][i+n]=1;
                a.m[i+n][i+n]=1;
            }

在草稿纸上画一下，这个会给你提示的。

我相信你已经在草稿纸上边滑了下，你会发现，这个在哪里出现了1，只有那个(2*n)*(2*n)的A矩阵里边，要怎么变成B，那么对应的右上，和右下正好都有一条线，就是左上右下对角线，对应的值就是1，其余的全部都是0，当然在你输入玩A之后会把左上边的那块全部赋值为对应的A矩阵的值，啊哈，大发现，这个不就是

|A    1|

|0    1|这个吗，太像了，不错，就是他，你太聪明了，么么哒！@——@

这个就是你要找到的B矩阵，这个时候只需要把B对应的矩阵来个n+1次方(这里为什么回事，不应该是n吗，嗯嗯，不过这里为了方便计算，就把他直接n+1，因为这样最后的结果就不需要再减去那个A^n，因为你这个时候得到的那个新的矩阵的右上方的那部分已经可以加到了A^n，而不再是A^n-1，懂了没？如果没动，可以QQ问我，我QQ2553627958，可以随时问我，学习就在于你自己主动不主动了)，得到的新的

(2*n)*(2*n)的矩阵，然后这个新的矩阵的右上方的值减个1就行，这就是最后的新的矩阵了也就是答案了，下边是我的代码：

#include<iostream>
#include<cstdio>
#include<cstring>
const int MAX=65;//这里为什么一定需要65？因为我直接开了(2*n)*(2*n)的矩阵，就是说

//对应的MAX的最小值就是n*2，也就是30*2=60.
using namespace std;
typedef struct
{
    int m[MAX][MAX];
}Matrax;
int mod,len;
Matrax multi(Matrax a,Matrax b)//矩阵a*b的算法
{
    Matrax c;
    int k,i,j;
    for(int i=1;i<=len;i++)
        for(int j=1;j<=len;j++)
    {
        c.m[i][j]=0;
        for(int k=1;k<=len;k++)
                c.m[i][j]=(c.m[i][j]+a.m[i][k]*b.m[k][j])%mod;//一定注意这里要取模，否则一定会WA
    }
    return c;
}
Matrax quick_pow(Matrax a,int n)//矩阵的奇数幂
{
    Matrax b;
    for(int i=1;i<=len;i++)//注意这里，这个就是矩阵的奇数幂的关键，和一般的奇数幂不同
        for(int j=1;j<=len;j++)
        b.m[i][j]=(i==j);
    while(n)//矩阵的快速幂
    {
        if(n&1)
            b=multi(b,a);
        n=n/2;
        a=multi(a,a);
    }
    return b;
}
int main()
{
    int n,k;
    while(~scanf("%d%d%d",&n,&k,&mod))
    {
        len=(n<<1);
        Matrax a;
        memset(a.m,0,sizeof a.m);
        for(int i=1;i<=n;i++)
            for(int j=1;j<=n;j++)
            scanf("%d",&a.m[i][j]);
            for(int i=1;i<=n;i++)
            {
                a.m[i][i+n]=1;
                a.m[i+n][i+n]=1;
            }
        a=quick_pow(a,k+1);
        for(int i=1;i<=n;i++)
            for(int j=n+1;j<=len;j++)
        {
            if(j==i+n)
                a.m[i][j]--;
//            a.m[i][i+n]--;//如果你把这个来替换掉上边的条件语句，那么你就和我翻过一样的而错误了，我找了好久，才发现这个错误，因为如果没有那个条件，那么也就代表这个语句每次都会执行一次，这样肯定错了，哎，心累啊。。。
            while(a.m[i][j]<0)//防止数据溢出
                a.m[i][j]+=mod;
        }
        for(int i=1;i<=n;i++)//最后打印出来
        {
            for(int j=n+1;j<len;j++)
            printf("%d ",a.m[i][j]);
            printf("%d\n",a.m[i][len]);
        }
    }
}