快速幂

快速幂是在进行底数相同的乘法运算而幂数又极大的情况下使用的一种算法。

模板

 https://www.luogu.com.cn/problem/P1226

实现

#include<iostream>
#include<cstdio>
#define int unsigned long long
using namespace std;
int a,b,q;
int qsm(int a,int b,int q)//快速幂
{
    int ans=1;
    while(b)
    {
        if(b&1)ans=(ans*a)%q;
        a=(a*a)%q;
        b>>=1;  //二进制运算，相当于b/=2;
    }
    return ans%q;
}
signed main()
{
    cin>>a>>b>>q;
    cout<<a<<"^"<<b<<" mod "<<q<<"="<<qsm(a,b,q);
    return 0;
}

分析原理，即a^b可以分解为a的2进制位数的次方，就是由X^m*Xⁿ = X^m+n

然后，当b&1（二进制运算，当第一位为1时返回true）时，说明这一位应该被ans乘，乘就对了。

如果要取模的话（对q取模），在每个数的后面取模就行了（防止数字过大）。

最后的最后，因为数字都比较大，记得开long long。

矩阵快速幂

这个算法的经典例题是求斐波那契数列第n项（mod q）的值。和快速幂一样，这类题指数非常的大。

要了解透彻矩阵快速幂，需要先了解一下矩阵乘法。这个知识在线性代数有所涉及（吧

若A为n×k矩阵，B为k×m矩阵，则它们的乘积AB(有时记做A·B)将是一个n×m矩阵。前一个矩阵的列数应该等于后一个矩阵的行数，得出的矩阵行数等于前一个矩阵的行数，列数等于后一个矩阵的行数。

举例：A、B均为3*3的矩阵：C=A*B，下面的代码会涉及到两种运算顺序，第一种就是直接一步到位求，第二种就是每次求一列，比如第一次，C00+=a00*b00,C01+=a00*b01……第二次C00+=a00*b10，C01+=a01*b11……以此类推。。。（网上摘录）

C00 = a00*b00 + a01*b10 + a02*b20
C01 = a00*b01 + a01*b11 + a02*b21 
C02 = a00*b02 + a01*b12 + a02*b22
C10 = a10*b00 + a11*b10 + a12*b20
C11 = a10*b00 + a11*b11 + a12*b21
C12 = a10*b02 + a11*b12 + a12*b22
C20 = a20*b00 + a21*b10 + a22*b20
C21 = a20*b01 + a21*b11 + a22*b21
C22 = a20*b02 + a21*b12 + a22*b22
C00 = a00*b00 + a01*b10 + a02*b20
C01 = a00*b01 + a01*b11 + a02*b21 
C02 = a00*b02 + a01*b12 + a02*b22

简而言之就是用前一个矩阵的行对应地乘下一个矩阵的列并把结果相加。用代码的方式可以表示为

for(int i=0;i<2;i++)  
        for(int j=0;j<2;j++)  
            for(int k=0;k<2;k++)  
        tmp.a[i][j]+=x.a[i][k]*y.a[k][j]　　,　　tmp.a[i][j]%=MOD;//如果题目让对某个数取模的话

这里的res矩阵用于临时存储结果，乘完后输出就行了。

矩阵乘法符合结合律但不符合交换律。

那么把矩阵当做运算的元素就有了矩阵快速幂，也就是乘幂次个矩阵。

例如斐波那契数列的问题。我们知道了矩阵乘法之后，可以发现

和快速幂的想法一样，用[1 1] [1 0]res这个矩阵一乘原来的矩阵i,i+1就可以得到i+1,i+2的矩阵了。然后和快速幂一样，只要每次让它乘自身就和快速幂一样了。

模板

https://www.luogu.com.cn/problem/P3390 //矩阵乘法

实现

 这个题是输入斐波那契数列的第1、2项，找第n项并对p取模。输入a1,a2,n,p。

matrix mul(matrix a,matrix b,ull q)//矩阵乘法
{
    matrix tmp;
    for(int i=0;i<=1;i++)
        for(int j=0;j<=1;j++)
        tmp.m[i][j]=0;
    for(int i=0;i<=1;i++)
        for(int j=0;j<=1;j++)
            for(int k=0;k<=1;k++)
            tmp.m[i][j]=(tmp.m[i][j]+(a.m[i][k]*b.m[k][j]%q)%q);
    return tmp;
}

ull fib(ull a1,ull a2,ull n,ull p)//解决问题
{
    ans.m[0][0]=a1,ans.m[0][1]=a2;//初始的两个数
    res.m[0][0]=1,res.m[0][1]=1,res.m[1][0]=1,res.m[1][1]=0;//初始化要乘的矩阵
    while(n)
    {
        if(n&1)ans=mul(ans,res,p);//如果要乘
        res=mul(res,res,p);
        n>>=1;
    }
    return ans.m[0][0];
}

int main()
{
    ull n,m,p,a1,a2;
    scanf("%lld%lld%lld%lld",&a1,&a2,&n,&p);//记得开longlong
    printf("%lld",fib(a1,a2,n-2,p));
    return 0;
}

至于矩阵乘法还能干什么，应该应用也挺多的，但是我还没那个水平……

大家一起进步吧XD

---

 

UPD：NOI2020出来了，发现D1T1用到了矩阵快速幂，于是来更新我写的水的一p的第一篇博客。

也就是更新一个例子：矩阵快速幂优化DP

我们先来个简单版的DP：给你一个图（的邻接矩阵），边权全为1，求给定时间从s到t的路径方案数。

记f[i]表示i时状态，f[i]=f[i-1]*g。

所以我们要求x直接矩阵快速幂乘g的幂次就行了。

 

---

 

2022.9.29

UPD：高中时不求甚解写的东西……也不是我不想了解，是竞赛要求我在理解之前就使用矩阵。突然感觉竞赛和文化课其实没有什么本质区别。

没时间大改，把错误的地方尽量删了。