题解 P3390 【【模板】矩阵快速幂】

· · 题解

矩阵快速幂这个东西嘛，其实很简单！

一句话:

矩阵快速幂=矩阵乘法+快速幂

所以，这篇题解就要从这两个东西说起：

1. 矩阵乘法

首先，不懂矩阵乘法的戳这里（如果你已经懂了，可以跳过）

还是一脸懵？？？我再来说一遍：

这回大家应该都懂了吧（如果图画的不好或者你还是没懂，请原谅）

那知道了怎么乘，到底怎么代码实现呢？？？

其实很简单：用模拟就可以啦！

大家先自己思考一下（一定要自己思考），思考完了再来看代码吧

5

4

3

2

1

相信大家代码都思考好了吧！接下来咱们来放一下代码：

    //本例中计算矩阵A×矩阵B，存到矩阵C里
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n;j++){  //找到c[i][j]，开始算
            for(int k=1;k<=n;k++){
                c[i][j]+=a[i][k]*b[k][j];//矩阵乘法的定义
            }
        }
    }

相信大家都已经get到矩阵乘法的原理和代码实现了吧！

顺便再提一句：主对角线（即左上-右下对角线）全是1，其它地方全是0的矩阵是单位矩阵也就是说单位矩阵长这样：
1 0 0
0 1 0
0 0 1

它就像数的乘法里1乘几就得几一样，它乘上哪个矩阵，就是哪个矩阵，为什么？自己乘一下看看就知道啦~

2.快速幂

首先，不会快速幂的看过来 （如果你AC了这道题或已学过相关内容，请跳过）

咱先想想这个问题：如果RMB是这样编排的：1元，2元，4元，8元......也就是面值都是2的几次方，那如果带齐“一套钱”，买几块的都能凑起来，不用找零了！为什么呢？？？

1=1 2=2 3=1+2 4=4 5=1+4 6=2+4 7=1+2+4 ......

真的哎！那为什么这么神奇呢？别急，既然是2的几次方，那咱要不....把每个数转化成二进制试试？？

1=1 2=10 3=11 4=100 5=101 6=110 7=111 ......

发现了吗？如果选了这个数，那对应的一位必然是1！好理解！二进制的第几位，表示的就是二的几次方嘛！而这些要么选（1），要么不选（0），必然有一种状态满足你！

我们的快速幂也是这个原理：对于任意一个数，都可以拆成2^k的数构成的和，所以求a^b，我们就可以把它拆成某几个a^{2^k}的积！

for example:

3^{11}=3^{8+2+1}=3^8\times3^2\times3^1

然后呢，把a每次都平方，a第一次平方之后会变成a^2，再平方就是a^4，......就可以快速访问每个a^{2^k}的状态了

处理的时候呢，把指数 b转化成二进制（就是除以2，膜2，再除以2，再膜2...），如果这一位是1，那就把当前的a乘到答案里去！

//求a的b次方，答案存到ans里
int ans=1;
while(b>0){
    if(b%2==1)  //如果当前位是1，那答案就包含这位，把a乘到ans里面
    ans=ans*a;
    a=a*a;//把a平方，切到下一状态
    b/=2;//除以2就不用说了吧
}

当然你可以用位运算优化，把/2改成>>1，把%2改成&1，据说会更快哦！

3.矩阵快速幂

如果你已了解并熟练掌握以上两点，那你离AC就只有一步之遥了

I have a 矩阵！

I have a 快速幂！

Ah~！矩阵快速幂！

没错，就这么简单！因为矩阵乘法是满足结合律的，所以一个单位矩阵用快速幂乘上n次就是我们的矩阵快速幂了，或者说，矩阵快速幂就是把快速幂里的乘法换成矩阵乘法

#include<cstdio>
using namespace std;
typedef long long ll; //十年OI一场空，不开long long见祖宗
const int syk=1e9+7;//syk是一个特别厉害的大佬
ll a[105][105],b;
int n;
ll ans[105][105]={0};
inline ll read(){
    char c=getchar();
    ll f=1,x=0;
    while(c<'0'||c>'9'){  //读入优化，这里就不解释了，想学的可以翻我之前的某篇题解
        if(c=='-')
            f=-1;
        c=getchar();
    }
    while(c>='0'&&c<='9'){
        x=(x<<1)+(x<<3)+(c^'0');
        c=getchar();
    }
    return x*f;
}
inline void jzcf1(){  //快速幂里的第一个乘法式子子ans=ans*a
    ll c[105][105]={0};
    for(int k=1;k<=n;k++){
        for(int i=1;i<=n;i++){
            for(int j=1;j<=n;j++){
                c[i][j]=(c[i][j]+ans[i][k]*a[k][j])%syk; //注意膜1e9+7
            }
        }
    }   
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n;j++){
            ans[i][j]=c[i][j];
        }
    }
}
inline void jzcf2(){//快速幂里的第而个乘法式子子a=a*a
    ll c[105][105]={0};
    for(int k=1;k<=n;k++){
        for(int i=1;i<=n;i++){
            for(int j=1;j<=n;j++){
                c[i][j]=(c[i][j]+a[i][k]*a[k][j])%syk;
            }
        }
    }   
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n;j++){
            a[i][j]=c[i][j];
        }
    }
}
int main(){
    n=(int)read(),b=read();
    for(register int i=1;i<=n;i++){
        for(register int j=1;j<=n;j++)
        a[i][j]=read();//读入
    }
    for(register int i=1;i<=n;i++)
    ans[i][i]=1;  //把ans初始化成单位矩阵
    while(b){  //快速幂来啦
        if(b&1)
        jzcf1();//把快速幂里的乘法改成矩阵乘法
        jzcf2();
        b>>=1;
    }
    for(register int i=1;i<=n;i++){ //输出
        for(register int j=1;j<=n;j++)
        printf("%lld ",ans[i][j]%syk);
        printf("\n");
    }
    return 0;
}

其实，只要掌握了矩阵乘法和快速幂的知识，矩阵快速幂其实一点都不难！屏幕前的小伙伴们，恭喜你们有get到了一个新的模板题！

The end ......