【数论】矩阵加速

刺骨的言语ヽ痛彻心扉 2021-12-14 14:47 239阅读 0赞

阅读本文之前 你应该有以下知识背景

[矩阵乘法法则][Link 1] && [矩阵快速幂][Link 2] && 小学数学水平 && 普及组代码能力 && 普及组思维能力

--------------------

矩阵加速的本质是一个用了玄学优化的递推

根据这个解释 我们可以首先知道它是一个递推 //23333

普通的递推是通过不断f\[i\] = f\[i-x\] (某种运算)f\[i-y\]来完成的

举个简单例子 斐波那契数列 f\[x\]=f\[x-1\]+f\[x-2\]

那么由于f\[x\]由f\[x-1\]+f\[x-2\]决定 所以我们在求f\[x\] 之前必须知道f\[x-1\]和f\[x-2\]

当然了普通人都是这么想的 但是到了矩阵手里这玩意就不是这样了

想像一个过程 一个矩阵乘了自己好几遍之后 再乘以一个已给出的矩阵 就是你的答案

对于斐波那契来说就是 (一个矩阵)^k \* (1,1) =answer

wow听起来很疯狂但是就是**矩阵加速**

--------------------

怎么做矩阵加速？

首先来看斐波那契的矩阵加速//注意 本人习惯使用竖着的f\[x\]\[2\]

你把递推边界放在一个矩阵里\{1,1\}然后乘以\{ \{0,1\}，\{1,1\}\}//这个怎么来的待会会说

得到\{1,2\}

你再试几组 就会发现每次都会往后走一个 即\{ f\[x-2\],f\[x-1\] \} --->\{ f\[x-1\],f\[x\] \}

你可能会说这不是比递推还慢？

但是递推能快速幂吗？ 抱歉 矩阵能

所以你要递推n次 就求一个\{ \{0,1\},\{1,1\}\}^n 然后把它乘以\{1,1\}

那么问题转化为 怎么构造转移矩阵

其实也就是在每个位置上体现递推规则

\[x-2\]->\[x-1\] 所以是\{0,1\}

\[x-1\]->\[x\]所以是\{1,1\}（\[x\] 和\[x+1\]相加）

OK讲解完毕 配几道题

斐波那契：

#include<bits/stdc++.h>
    #define ll long long
    #define p 1000000007
    
    using namespace std;
    ll k,n,s;
    struct mat{
        ll m[3][2];
    }f;//ans列 
    struct smat{
        ll m[3][3];
    }spd,e;
    
    smat smul(smat x,smat y){
        smat r;
        for(int i=1;i<=2;i++){
            for(int j=1;j<=2;j++){
                r.m[i][j] = 0;
            }
        }
        for(int i=1;i<=2;i++){
            for(int j=1;j<=2;j++){
                for(int k=1;k<=2;k++){
                    r.m[i][j] += x.m[i][k] * y.m[k][j] %p;
                    r.m[i][j] %=p;
                }
            }
        }
        return r;
    }
    
    mat mul (smat x,mat y){
        mat r;
        for(int i=1;i<=2;i++){
            r.m[i][1]=0;
        }
        for(int i=1;i<=2;i++){
            for(int k=1;k<=2;k++){
                r.m[i][1] += x.m[i][k] * y.m[k][1] %p;
                r.m[i][1] %=p;
            }
        }
        return r;
    }
    
    smat pow(smat x,ll k){
        smat ans = e;
        while (k>0){
            if(k&1){
                ans=smul(ans,x);
            }
            x=smul(x,x);
            k>>=1;
        }
        return ans;
    }
    
    int main(){
        spd.m[1][2] = spd.m[2][1] = spd.m[2][2] = 1;
        e.m[1][1] = e.m[2][2] = 1;
        f.m[2][1] = f.m[1][1] = 1;
        cin>>n;
        smat vf=pow(spd,n-2);
        f=mul(vf,f);
        cout<<f.m[2][1];
        return 0;
    }

板子[ P1939 【模板】矩阵加速（数列）][P1939]

其实大同小异吧。。。自己看代码

#include<bits/stdc++.h>
    #define ll long long
    #define p 1000000007
    
    using namespace std;
    ll k,n,s;
    struct mat{
        ll m[4][2];
    }f;//ans列 
    struct smat{
        ll m[4][4];
    }spd,e;
    
    smat smul(smat x,smat y){
        smat r;
        for(int i=1;i<=3;i++){
            for(int j=1;j<=3;j++){
                r.m[i][j] = 0;
            }
        }
        for(int i=1;i<=3;i++){
            for(int j=1;j<=3;j++){
                for(int k=1;k<=3;k++){
                    r.m[i][j] += x.m[i][k] * y.m[k][j] %p;
                    r.m[i][j] %=p;
                }
            }
        }
        return r;
    }
    
    mat mul (smat x,mat y){
        mat r;
        for(int i=1;i<=3;i++){
            r.m[i][1]=0;
        }
        for(int i=1;i<=3;i++){
            for(int k=1;k<=3;k++){
                r.m[i][1] += x.m[i][k] * y.m[k][1] %p;
                r.m[i][1] %=p;
            }
        }
        return r;
    }
    
    smat pow(smat x,ll k){
        smat ans = e;
        while (k>0){
            if(k&1){
                ans=smul(ans,x);
            }
            x=smul(x,x);
            k>>=1;
        }
        return ans;
    }
    
    int main(){
        spd.m[1][2] = spd.m[2][3] = spd.m[3][1] = spd.m[3][3] = 1;
        e.m[1][1] = e.m[2][2] = e.m[3][3] = 1;
        f.m[3][1] = f.m[2][1] = f.m[1][1] = 1;
        ll wow;
        cin>>wow;
        for(int i=1;i<=wow;i++){
            cin>>n;
            
            smat vf=pow(spd,n-3);
            f=mul(vf,f);
            cout<<f.m[3][1]<<endl;
            f.m[3][1] = f.m[2][1] = f.m[1][1] = 1;
        }
        return 0;
    }

本文是针对有快速幂基础的

所以没有那么多解释

# TAG：SIN\_XIII ⑨ #

转载于:https://www.cnblogs.com/SINXIII/p/10543427.html

[Link 1]: https://baike.baidu.com/item/%E7%9F%A9%E9%98%B5%E4%B9%98%E6%B3%95/5446029?fr=aladdin
[Link 2]: https://www.luogu.org/problemnew/show/P3390
[P1939]: https://www.luogu.org/problemnew/show/P1939