矩阵乘法

矩阵

矩阵的定义:

一个n*m的矩阵可以看作是一个二维数组

设A是\(n * m\)矩阵,B是 \(m * p\)矩阵

则C就是\(n * p\) 矩阵 并且

\[\sum _{k=1}^{m} A_{i,k}*B_{k,j} \]

矩阵乘法满足结合律,即\((A*B)*C\) =\(A*(B*C)\)
满足分配律,即\((A+B)*C\)=\(A*C+B*C\)


矩阵的递推应用

若F1 是 \(1*n\)矩阵,F2是\(n*n\)矩阵,F3=F1*F2

F1,F3可看作一维数组

这时我们发现,矩阵乘法就可以用来做一些递推式

举个image

假如一个递推式是 \(f[i]=f[i-1]+f[i-2]\)(斐波那契数列)

那么它的我们就可以设它的A矩阵为

A={\(f[i-1]\) \(f[i]\)}

想要得到的C矩阵就应该是

C={\(f[i]\) \(f[i+1]\)}

此时,我们是不是只要找到一个\(n*n\)的B矩阵就可以使得\(A*B=C\)

那么怎么构造B矩阵呢?首先我们要理解一下C数组

\(C[i][j]\)表示的应该就是A数组第i行的每个数分别乘上B数组第j行的每个数

那么对于这个题来说

我们要得到的\(C[1][1]\)是不是\(f[i]\)

那么我们只要\(f[i-1]\)乘上0,f[i]乘上1是不是就行了

所以说B的第一列是0,1

另外

我们观察\(f[i+1]\) 显然\(f[i+1]=f[i]+f[i-1]\)

所以我们要让\(f[i-1]*1+f[i]*1\)

所以B数组就是

0 1

1 1

我们可以这样写
同时我们发现,F3=F1*F2

F4=F3*F2

\(F5=F4*F2=F3*F2*F2=F1*F2^3\)

\(Fn=F1*F2^{n-1}\)

upd:我发现了一种新的推式子的好方法:
对于

\[a_{x}= \begin{cases}1 & x \in\{1,2,3\} \\ a_{x-1}+a_{x-3} & x \geq 4\end{cases} \]

我们可以选择进行一次新的矩阵变换:
\(A={{f_{i-3},f_{i-2},f_{i-1}}}\)
\(C={{f_{i-2},f_{i-1},f_{i}}}\)
我们尝试竖着写(等会会很有用的)

\[C= \left[\begin{array}{c} f_i \\ f_{i-1} \\ f_{i-2} \end{array}\right] \]

\[A= \left[\begin{array}{c} f_{i-1} \\ f_{i-2} \\ f_{i-3} \end{array}\right] \]

那么我们这样做:
直接看看C的每一个式子是怎么推出来的,但是位置和顺序不能换,0也得写上去

\[f_{i}=1\times f_{i-1}+0\times f_{i-2}+1\times f_{i-3};\\ f_{i-1}=1\times f_{i-1}+0\times f_{i-2}+0\times f_{i-3};\\ f_{i-2}=0\times f_{i-1}+1\times f_{i-2}+0\times f_{i-3}; \]

把前面的每一个系数取出来
所以我们的矩阵 \(B\) 最后就长成

\[1~~0~~1\\ 1~~0~~0\\ 0~~1~~0 \]

这个是对的,不信自己试试看?

但是应该你得把A和B换一下,要不不满足矩阵乘了

最后贴一份矩阵快速幂的代码罢

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#include<queue>
#include<cctype>
#include<bitset>
#include<vector>
#include<map>
#include<set>
#define int long long
using namespace std;
const int N=110;
const int mod=1e9+7;
typedef pair<int,int> P;
namespace scan{
    template <typename T>
    void read(T &x){
        x=0; char c=getchar(); int f=0;
        for(; !isdigit(c); c=getchar()) f |= c=='-';
        for(; isdigit(c); c=getchar()) x = x*10+(c^48);
        if(f) x = -x;
    }

    template <typename T>
    void write(T x,char ch){
        if(x<0) putchar('-'), x = -x;
        static short st[30],tp;
        do st[++tp] = x%10, x/=10; while(x);
        while(tp) putchar(st[tp--] | 48);
        putchar(ch);
    }
}
using namespace scan;
int n,k;
struct Matrix{
    int a[N][N];
    Matrix(){memset(a,0,sizeof(a));}
    inline void build(){
        for(int i=1;i<=n;i++) a[i][i]=1;
    }
}A;

Matrix operator * (const Matrix &A,const Matrix &B){
    Matrix C;
    for(int k=1;k<=n;k++)
        for(int i=1;i<=n;i++)
            for(int j=1;j<=n;j++)
                C.a[i][j]=(C.a[i][j]+A.a[i][k]*B.a[k][j])%mod;
    return C;
}
Matrix ans;
inline void qpow(int k){
    while(k){
        if(k&1) ans=ans*A;
        A=A*A;
        k>>=1;
    }
}
signed main(){
    read(n),read(k);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++) read(A.a[i][j]);
    ans.build();
    qpow(k);
    for(int i=1;i<=n;i++){
        for(int j=1;j<=n;j++)
            write(ans.a[i][j],' ');
        putchar('\n');   
    }
    return 0;
}
posted @ 2021-09-23 22:00  RevolutionBP  阅读(238)  评论(0编辑  收藏  举报