模板矩阵加速

Posted winlere

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了模板矩阵加速相关的知识,希望对你有一定的参考价值。

矩阵加速,专门用来解决一些递推的关系,其原理和矩阵运算的法则有关

由于矩阵的乘法有交换律和结合律,所以我们可以通过矩阵快速幂来快速求解递推关系,一般时间复杂度是O(nlogn)。

矩阵快速幂很简单,写一下模板就会了,但是推导单位矩阵是个难题。

一般地,我们推导单位矩阵时,有这几个步骤。

  1. 确定递推初始条件
  2. 确定递推式子的系数
  3. 确定单位矩阵的大小,一般来说和递推式涉及到的已知量个数有关,必须考虑系数为0的原量,不能忽略
  4. 填单位矩阵,注意一下几点:
    1. 第一行才是生成新数据的力量。
    2. 所以后面的行只需要把自己表示出来就好了。(为什么是用“行”来表示已知量?根据矩阵运算法则)

一般容易犯的错误:

  • 递推关系写错了
  • 单位矩阵错误了
  • 矩阵乘法写错了
  • 答案输出错了,没有找到到底应该输出哪一个
  • 矩阵的指数确定错了
  • %错了
  • 没开long long导致乘法爆了

下贴斐波那契数列矩阵加速版本

#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
#define RP(t,a,b) for(ll t=(a),edd=(b);t<=edd;t++)
const int mod=1e9+7;
typedef unsigned long long ll;
const int n=2;
ll k;
struct mtx {
    ll data[3][3];
    ll* operator [](int x){return data[x];}
    mtx() {
        memset(data,0,sizeof data);
    }
    mtx operator *(const mtx& x) {
        mtx temp;
        RP(t,0,n-1)
        RP(i,0,n-1)
        RP(k,0,n-1)
        temp.data[t][i]=(temp.data[t][i]+(data[t][k]*x.data[k][i])%mod)%mod;
        return temp;
    }
    mtx operator *=(const mtx& x) {
        return (*this)=(*this)*x;
    }
    mtx operator ^(const ll& p) {
        ll b=p;
        mtx ans,base;
        ans.unis();
        base=(*this);
        while(b) {
            if(b&1)
                ans*=base;
            base*=base;
            b>>=1;
        }
        return ans;
    }
    mtx operator ^=(const ll& p) {
        return (*this)=(*this)^p;
    }
    void unis() {
        for(ll t=0;t<n;t++)
            data[t][t]=1;
    }
    void print() {
        for(ll t=0; t<n; t++) {
            for(ll ti=0; ti<n; ti++)
                cout<<data[t][ti]<<‘ ‘;
            cout<<"
";
        }
        return;
    }
};
int T;
int main() {
    //      freopen("in.in","r",stdin);
    cin>>k;
    mtx qaq;
    qaq.data[0][1]=qaq.data[1][0]=qaq.data[1][1]=1;
    qaq^=k-1;
    ll asa=qaq.data[0][0]%mod+qaq.data[1][0]%mod;
    cout<<asa%mod<<endl;
    return 0;
}

然后又附luogu矩阵加速板子的代码

#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
#define RP(t,a,b) for(ll t=(a),edd=(b);t<=edd;t++)
const int mod=1e9+7;
typedef long long ll;
const int n=3;
int k;
struct mtx {
    ll data[3][3];
    mtx() {
        memset(data,0,sizeof data);
    }
    mtx operator *(const mtx& x) {
        mtx temp;
        RP(t,0,n-1)
        RP(i,0,n-1)
        RP(k,0,n-1)
        temp.data[t][i]=(temp.data[t][i]+(data[t][k]*x.data[k][i])%mod)%mod;
        return temp;
    }
    mtx operator *=(const mtx& x) {
        return (*this)=(*this)*x;
    }
    mtx operator ^(const ll& p) {
        ll b=p;
        mtx ans,base;
        ans.unis();
        base=(*this);
        while(b) {
            if(b&1)
                ans*=base;
            base*=base;
            b>>=1;
        }
        return ans;
    }
    mtx operator ^=(const ll& p) {
        return (*this)=(*this)^p;
    }
    void unis() {
        for(ll t=0;t<n;t++)
            data[t][t]=1;
    }
    void print() {
        for(ll t=0; t<n; t++) {
            for(ll ti=0; ti<n; ti++)
                cout<<data[t][ti]<<‘ ‘;
            cout<<"
";
        }
        return;
    }
};
int T;
int main() {
    //  freopen("in.in","r",stdin);
    cin>>T;
    while(T--){
        cin>>k;
        mtx qaq;
        mtx temp;
        qaq.data[0][0]=qaq.data[0][2]=qaq.data[1][0]=qaq.data[2][1]=1;
        qaq^=k;
        int pos=0;
        cout<<qaq.data[1][0]<<endl;
    }
    return 0;
}
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
#define RP(t,a,b) for(ll t=(a),edd=(b);t<=edd;t++)
const int mod=1e9+7;
typedef long long ll;
const int n=3;
int k;
struct mtx {
    ll data[3][3];
    mtx() {
        memset(data,0,sizeof data);
    }
    mtx operator *(const mtx& x) {
        mtx temp;
        RP(t,0,n-1)
        RP(i,0,n-1)
        RP(k,0,n-1)
        temp.data[t][i]=(temp.data[t][i]+(data[t][k]*x.data[k][i])%mod)%mod;
        return temp;
    }
    mtx operator *=(const mtx& x) {
        return (*this)=(*this)*x;
    }
    mtx operator ^(const ll& p) {
        ll b=p;
        mtx ans,base;
        ans.unis();
        base=(*this);
        while(b) {
            if(b&1)
                ans*=base;
            base*=base;
            b>>=1;
        }
        return ans;
    }
    mtx operator ^=(const ll& p) {
        return (*this)=(*this)^p;
    }
    void unis() {
        for(ll t=0;t<n;t++)
            data[t][t]=1;
    }
    void print() {
        for(ll t=0; t<n; t++) {
            for(ll ti=0; ti<n; ti++)
                cout<<data[t][ti]<<‘ ‘;
            cout<<"
";
        }
        return;
    }
};
int T;
int main() {
    //  freopen("in.in","r",stdin);
    cin>>T;
    while(T--){
        cin>>k;
        mtx qaq;
        mtx temp;
        qaq.data[0][0]=qaq.data[0][2]=qaq.data[1][0]=qaq.data[2][1]=1;
        qaq^=k;
        int pos=0;
        cout<<qaq.data[1][0]<<endl;
    }
    return 0;
}

唉我太弱了orz

以上是关于模板矩阵加速的主要内容,如果未能解决你的问题,请参考以下文章

Luogu P3390 模板矩阵快速幂&&P1939 模板矩阵加速(数列)

P1939模板矩阵加速(数列)

P1939 模板矩阵加速(数列)

模板矩阵加速(数列)

P1939 模板矩阵加速(数列)

模板矩阵加速(数列)