矩阵小专题(矩阵快速幂+矩阵加速)

1.什么是矩阵?

矩阵(数学术语)_百度百科


2.矩阵快速幂

首先要知道,只有n*n的矩阵能乘以自身(否则不符合矩阵相乘的条件)

然后要明白普通的快速幂的原理(本质是把幂次二分,代码如下)

inline int Pow(int a,int b){
    long long res=1, bs=a;
    while(b>0){if(b&1) res=res*bs%mod; bs=bs*bs%mod; b>>=1;}
    return res;
}

这里最终答案res从1开始,逐步让res乘以base值,以及让base值翻倍,从而达到了二分幂次的效果。

换成矩阵,思路完全一样。但是要进行两个替换:

  1. 把 res=1 换成 res等于一个单位矩阵
  2. 数字相乘换成矩阵相乘,这一步可以用运算符重载来完成。

 附上模板题链接和代码:【模板】矩阵快速幂 - 洛谷

#include <bits/stdc++.h>
using namespace std;
#define FOR(i, a, b) for (int i = (a); i <= (b); i++)
#define int long long
#define pii pair<int,int>
const int N = 105, mod=1e9+7;
int n, k; //矩阵大小,矩阵次方
struct juzhen{
    int a[N][N];
    juzhen(){memset(a,0,sizeof(a));} //初始为0矩阵
    inline void build(){FOR(i,1,n) a[i][i]=1;} //建立单位矩阵
    void print(){
        FOR(i,1,n){
            FOR(j,1,n) cout<<a[i][j]<<' ';
            cout<<'\n';
        }
    }
    juzhen operator*(const juzhen& b){ //重载矩阵乘法
        juzhen res;
        FOR(i,1,n) FOR(j,1,n){
            FOR(k,1,n) res.a[i][j] = (res.a[i][j] + this->a[i][k]*b.a[k][j]%mod)%mod;
        }
        return res;
    }
}a;

void solve(){
    cin>>n>>k;
    FOR(i,1,n) FOR(j,1,n) cin>>a.a[i][j];
    juzhen ans; ans.build(); //从一个单位矩阵开始
    while(k>0){ //矩阵快速幂计算
        if(k&1) ans=ans*a;
        a = a*a; k>>=1;
    }
    ans.print();
}
signed main(){
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    int T=1; //cin>>T;
    while(T--) solve();
}

3.矩阵加速 

矩阵加速可以用来计算一些dp的计算,加速的原理就是前面说到的矩阵快速幂。

1) 首先我们以一个简单的问题引入:斐波那契数列 - 洛谷

这就是普通的斐波那契数列,但是n值达到了恐怖的2e9,显然原本O(n)的普通dp思路不可行。

那么,我们就要用到矩阵加速了。

考虑构造一个1*2的矩阵[f[i-1],f[i]],用这个矩阵进行递推。

由斐波那契的递推式f[i] = f[i-2]+f[i-1]

得到矩阵的递推式[f[i-1],f[i]] = [f[i-2],f[i-1]]*\begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix}

再得到最终的计算式[f[n-1],f[n]] = [f[0],f[1]]*\begin{bmatrix} 0 & 1\\ 1 & 1 \end{bmatrix} ^{n-1}

最终代码如下:

#include <bits/stdc++.h>
using namespace std;
#define FOR(i, a, b) for (int i = (a); i <= (b); i++)
#define int long long
#define pii pair<int,int>
const int mod = 1e9+7;
int n;

struct node{
    int a[3][3];
    node(){memset(a,0,sizeof(a));} //初始为0矩阵
    void build(){a[1][1]=0,a[1][2]=a[2][1]=a[2][2]=1;} //建立变换矩阵
    node operator*(const node&b){
        node res;
        FOR(i,1,2) FOR(j,1,2){
            FOR(k,1,2) res.a[i][j] = (res.a[i][j]+a[i][k]*b.a[k][j]%mod)%mod;
        }
        return res;
    }
};
node Pow(node a, int b){
    node res; res.a[1][1]=res.a[2][2]=1; //建立res为单位矩阵
    while(b>0){ //矩阵快速幂
        if(b&1) res=res*a;
        a=a*a; b>>=1;
    }
    return res;
}
pii cheng(pii a, node b){
    pii res;
    res.first = (a.first*b.a[1][1]%mod + a.second*b.a[2][1]%mod)%mod;
    res.second = (a.first*b.a[1][2]%mod + a.second*b.a[2][2]%mod)%mod;
    return res;
}
void solve(){
    cin>>n;
    node a; a.build();
    a = Pow(a,n-1);
    cout<<cheng({0,1},a).second<<'\n';
}
signed main(){
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    int T=1; //cin>>T;
    while(T--) solve();
}

2) 再推广到一般的矩阵加速,也是一样的做法,看下这道题:【模板】矩阵加速(数列) - 洛谷

容易得到矩阵转移式为:[a,b,c]\rightarrow [b,c,a+c] 

所以这次的转换矩阵就变成了\begin{bmatrix} 0 & 0 & 1\\ 1 & 0 & 0\\ 0 & 1 & 1 \end{bmatrix}

最终代码如下:

#include <bits/stdc++.h>
using namespace std;
#define FOR(i, a, b) for (int i = (a); i <= (b); i++)
#define int long long
#define pii pair<int,int>
const int mod = 1e9+7;
int n;

struct node{
    int a[4][4];
    node(){memset(a,0,sizeof(a));} //初始为0矩阵
    void build(){a[1][3]=a[2][1]=a[3][2]=a[3][3]=1;} //建立变换矩阵
    node operator*(const node&b){
        node res;
        FOR(i,1,3) FOR(j,1,3){
            FOR(k,1,3) res.a[i][j] = (res.a[i][j]+a[i][k]*b.a[k][j]%mod)%mod;
        }
        return res;
    }
    void print(){
        FOR(i,1,3){
            FOR(j,1,3) cout<<a[i][j]<<' ';
            cout<<'\n';
        }
    }
};
node Pow(node a, int b){
    node res; res.a[1][1]=res.a[2][2]=res.a[3][3]=1; //建立res为单位矩阵
    while(b>0){ //矩阵快速幂
        if(b&1) res=res*a;
        a=a*a; b>>=1;
    }
    return res;
}
struct tii{
    int a[4];
    tii(){memset(a,0,sizeof(a));}
    tii(int x,int y,int z){a[1]=x,a[2]=y,a[3]=z;}
};
tii cheng(tii a,node b){
    tii res;
    FOR(i,1,3){
        FOR(j,1,3) res.a[i] = (res.a[i] + a.a[j]*b.a[j][i]%mod)%mod;
    }
    return res;
}
void solve(){
    cin>>n;
    node a; a.build();
    a = Pow(a,n-1);
    cout<<cheng(tii(0,0,1),a).a[3]<<'\n';
}
signed main(){
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    int T=1; cin>>T;
    while(T--) solve();
}

  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值