码迷,mamicode.com
首页 > 其他好文 > 详细

【模板】矩阵加速

时间:2019-01-23 12:28:54      阅读:182      评论:0      收藏:0      [点我收藏+]

标签:end   ons   const   typedef   cstring   ring   注意   opera   main   

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

由于矩阵的乘法有交换律和结合律,所以我们可以通过矩阵快速幂来快速求解递推关系,一般时间复杂度是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<<"\n";
        }
        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<<"\n";
        }
        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<<"\n";
        }
        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

【模板】矩阵加速

标签:end   ons   const   typedef   cstring   ring   注意   opera   main   

原文地址:https://www.cnblogs.com/winlere/p/10308119.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!