【模板】【代数】矩阵乘法和矩阵快速幂

矩阵乘法

定义

设A为m*p的矩阵，B为p*n的矩阵，那么称m*n的矩阵C为矩阵A与B的乘积，记作$C=AB$，其中矩阵C中的第 i 行第 j 列元素可以表示为：

$$(AB)_{ij}=\sum_{k=1}^pa_{ik}b_{kj}=a_{i1}b_{1j}+a_{i2}b_{2j}+\cdots+a_{ip}b_{pj}$$

注意事项

1.当矩阵 A 的列数等于矩阵 B 的行数时， A 与 B 可以相乘。
2.矩阵 C 的行数等于矩阵 A 的行数， C 的列数等于 B 的列数。
3.乘积 C 的第 m 行第 n 列的元素等于矩阵 A 的第 m 行的元素与矩阵 B 的第 n 列对应元素乘积之和。

基本性质

1.乘法结合律： $(AB)C=A(BC)$
2.乘法左分配律：$(A+B)C=AC+BC$
3.乘法右分配律：$C(A+B)=CA+CB$
4.对数乘的结合性$k(AB)=(kA)B=A(kB)$
5.矩阵乘法一般不满足交换律

代码

struct matrix{
    LL z[mxl][mxr];
    int m,n;//矩阵的行数和列数 
};
matrix mul(matrix x,matrix y)//矩阵乘法 
{
    matrix t;
    t.m=x.m;//新矩阵的行数等于第一个矩阵的行数 
    t.n=y.n;//新矩阵的列数等于第二个矩阵的列数 
    memset(t.z,0,sizeof(t.z));
    for(int i=1;i<=x.m;i++)
        for(int j=1;j<=y.n;j++)
            for(int k=1;k<=y.m;k++)
                t.z[i][j]+=x.z[i][k]*y.z[k][j];
    return t;
}

时间复杂度$O(n^3)$

矩阵快速幂

与数的快速幂原理相同，只不过运算对象换成了矩阵。

数的快速幂请看http://blog.csdn.net/george__yu/article/details/77259517

矩阵快速幂能在$O(n^3log_2k)$的时间内算出$A^k$，其中A是一个矩阵。

代码

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long LL;
const int mxl=103,mxr=103;
struct matrix{
    LL z[mxl][mxr];
    int m,n;//矩阵的行数和列数 
}res,ori;//res:单位矩阵 ,ori:输入的矩阵 
LL p=1000000007,k;
void init()
{
    scanf("%d%d%lld",&ori.m,&ori.n,&k);//输入指数，行数，列数 
    res.m=ori.m;res.n=ori.n;//初始化单位矩阵的行数和列数 
    for(int i=0;i<=res.m+1;i++) res.z[i][i]=1;//初始化单位矩阵 
    for(int i=1;i<=ori.m;i++)//输入矩阵 
        for(int j=1;j<=ori.n;j++)
            scanf("%lld",&ori.z[i][j]);
}
matrix mul(matrix x,matrix y)//矩阵乘法 
{
    matrix t;
    t.m=x.m;//新矩阵的行数等于第一个矩阵的行数 
    t.n=y.n;//新矩阵的列数等于第二个矩阵的列数 
    memset(t.z,0,sizeof(t.z));
    for(int i=1;i<=x.m;i++)
        for(int j=1;j<=x.n;j++)
        {
            t.z[i][j]=0;
            for(int k=1;k<=y.n;k++)
                t.z[i][j]=(t.z[i][j]+x.z[i][k]*y.z[k][j])%p;
        }
    return t;
}
void MatrixModpow(LL k)//矩阵快速幂 
{
    while(k)
    {
        if(k&1) res=mul(res,ori);
        ori=mul(ori,ori);
        k>>=1;
    }
}
void print(matrix x)//打印矩阵 
{
    for(int i=1;i<=x.m;i++)
    {
        for(int j=1;j<=x.n;j++)
            printf("%lld ",x.z[i][j]);
        printf("\n");
    }
}
int main()
{
    init();
    MatrixModpow(k);
    print(res);
    return 0;
}

再补一个代码

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cstdlib>
#include <cmath>

using namespace std;
typedef long long LL;
const int mxl = 103,mxr = 103;
const LL mod = 1e9+7;

struct Matrix{
    LL l,r;
    LL num[mxl][mxr];
    Matrix() {memset(num,0,sizeof(num)); l=0; r=0;}

    Matrix operator * (const Matrix &b) const
    {
        if(r!=b.l) return *this;
        Matrix c;
        c.l=l;c.r=b.l;
        for(int i=1;i<=l;i++)
            for(int j=1;j<=b.r;j++)
                for(int k=1;k<=r;k++)
                    c.num[i][j]=(c.num[i][j]+(num[i][k]*b.num[k][j])%mod)%mod;
        return c;   
    }

    Matrix ModPow(LL b)
    {
        Matrix res,a;
        a=*this;
        res.l=res.r=a.l;    
        for(int i=1;i<=res.l;i++) res.num[i][i]=1;
        while(b)
        {
            if(b&1) res=res*a;
            a=a*a;
            b>>=1;
        }
        return res;
    }

    inline void read(int ll,int rr)
    {
        l=ll;r=rr;
        for(int i=1;i<=l;i++)
            for(int j=1;j<=r;j++)
                scanf("%lld",&num[i][j]);
    }

    inline void write()
    {
        for(int i=1;i<=l;i++)
        {
            for(int j=1;j<=r;j++) printf("%lld ",num[i][j]);
            printf("\n");
        }
    }
};

LL n,k;
Matrix a,b;

int main()
{
    scanf("%lld%lld",&n,&k);
    a.read(n,n);
    a.ModPow(k).write();
    return 0;
}

用途：优化一类线性递推问题。

具体请看我的这一篇博客：http://blog.csdn.net/george__yu/article/details/77249237