矩阵快速幂

· · 个人记录

矩阵快速幂，就是矩阵乘法和快速幂的结合，所以在讲矩阵快速幂之前，先讲讲矩阵乘法和快速幂。

1. 矩阵

1.1 矩阵的定义

什么是矩阵？由m \times n个数a_{i,j}(i=1,2,\cdots,m;j=1,2,\cdots,n)排成的一个m行n列的数表

\begin{matrix} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n}\end{matrix}

被称为m行n列矩阵，简称m \times n矩阵。为了表示这是个整体总是加一个括弧，并用大写字母表示她，计做：

A = \left(\begin{matrix} a_{1,1} & a_{1,2} & \cdots & a_{1,n} \\ a_{2,1} & a_{2,2} & \cdots & a_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1} & a_{m,2} & \cdots & a_{m,n}\end{matrix}\right)

这个矩阵中的m \times n个数被称为矩阵A的元素，简称元。矩阵A也可被计做A_{m \times n}。

行矩阵（行向量）： 只有一行的矩阵。
列矩阵（列向量）： 只有一列的矩阵。
同型矩阵： 两个矩阵A、B的行数和列数都相等
单位矩阵： 在矩阵乘法中起着很大的作用，相当于乘法中的1，她是个方阵（即行数和列数相同的矩阵），左上角到右下角的对角线（被称为主对角线）上的元素都是1，除此之外都是0。如下就是一个4 \times 4的单位矩阵：

\left(\begin{matrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1\end{matrix}\right)

1.2 矩阵加法

我们定义两个m \times n的矩阵A和B相加为：

A+B = \left(\begin{matrix} a_{1,1}+b_{1,1} & a_{1,2}+b_{1,2} & \cdots & a_{1,n}+b_{1,n} \\ a_{2,1}+b_{2,1} & a_{2,2}+b_{2,2} & \cdots & a_{2,n}+b_{2,n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m,1}+b_{m,1} & a_{m,2}+b_{m,2} & \cdots & a_{m,n}+b_{m,n}\end{matrix}\right)

注意：

1.3 矩阵乘法

我们定义两个矩阵A_{m \times s}和B_{s \times n}相乘得到一个矩阵C_{m,n}，并且

c_{i,j} = a_{i,1}b_{1,j}+a_{i,2}b_{2,j} \cdots +a_{i,s}b_{s,j}=\sum^{s}_{k=1}a_{i,k}b_{k,j}

并将此计做C=AB。
注意： 矩阵乘法不满足交换律。

\tiny{\text{这里没放图是因为放了也看不懂（逃}}

如果还不理解的话可以看这里
用代码写出来就是下面这个亚子：

struct matrix
{
    int n,m;//n是行数，m是列数
    ll e[105][105];
}a;
matrix mul(matrix a,matrix b)
{
    matrix ans;
    ans.n=a.n;
    ans.m=b.m;
    for(int i=1;i<=ans.n;i++)
        for(int j=1;j<=ans.m;j++)
            ans.e[i][j]=0;//初始化
    for(int i=1;i<=a.n;i++)
        for(int j=1;j<=b.m;j++)
            for(int k=1;k<=a.m;k++)
                ans.e[i][j]=(ans.e[i][j]+(a.e[i][k]*b.e[k][j])%mod)%mod;//矩阵乘法结果可能很大，所以要取模
    return ans;
}

2. 快速幂

一般我们求x^y的方法是一个个乘过去，但这种方法太慢了。我们需要一种新的，快速的方法来求，而这种方法就是快速幂。这里我们举例来理解，就比如求3^{5}，用ans记录答案。5转成二进制就是{(101)}_2，那么我们从右往左开始枚举，首先，有个1，那么我们就记录ans=ans \times 3=1 \times 3=3，然后让3变成自己的平方，也就是3^2=9，然后接下来是个0，那就不用乘，然后再让9变成自己的平方81，然后又是个1，所以ans=ans \times 81=3 \times 81=243。最后，得到3^5 = 243。手算一下可以发现这是对的。当然x^y很可能是个很大的数，所以一般都会取模。
代码实现长这个亚子：

long long power(long a,long b,long c)
{
    long long ans=1;
    while(b)//没有枚举完
    {
        if(b&1) ans=(ans*a)%c;//这一位是否是1
        a=(a*a)%c;
        b>>=1;//走到下一位
    }
    return ans%c;
}

3. 矩阵快速幂

讲完了上面这些，那矩阵快速幂就很好理解了！就是把运算改成矩阵运算。
首先，ans和返回值都应该是一个矩阵；然后，ans应该初始化为单位矩阵，原因上面已经讲了；（1.1 矩阵的定义）再然后，里面的乘法应该改成矩阵乘法。最后提醒一句，矩阵快速幂必须是个方阵，否则要两个同型矩阵（自己乘自己，当然是同型矩阵啦~）可以实现矩阵乘法是不可能的。
代码实现：

struct matrix
{
    int n,m;
    ll e[105][105];
}a;
matrix mul(matrix a,matrix b)
{
    matrix ans;
    ans.n=a.n;
    ans.m=b.m;
    for(int i=1;i<=ans.n;i++)
        for(int j=1;j<=ans.m;j++)
            ans.e[i][j]=0;
    for(int i=1;i<=a.n;i++)
        for(int j=1;j<=b.m;j++)
            for(int k=1;k<=a.m;k++)
                ans.e[i][j]=(ans.e[i][j]+(a.e[i][k]*b.e[k][j])%mod)%mod;//快速幂的取模就在这里了！
    return ans;
}//代码和上面的区别
matrix power(matrix a,ll b)
{
    matrix ans;
    ans.n=a.n;
    ans.m=a.m;
    for(int i=1;i<=ans.n;i++)
        for(int j=1;j<=ans.m;j++)
            if(i==j) ans.e[i][j]=1;
            else ans.e[i][j]=0;//初始化为单元矩阵
    while(b)
    {
        if(b&1) ans=mul(ans,a);//要把乘换成矩阵乘法哦！
        a=mul(a,a);//这里也一样呢
        b>>=1;
    }//快速幂板子
    return ans;
}

顺便挂上矩阵加速的友链