矩陣快速冪

矩陣優化能夠常常利用在遞推式中。ios

首先了解一下矩陣乘法的法則。優化

\(\begin{bmatrix}a&b\\c&d\end{bmatrix}\) \(\times\) \(\begin{bmatrix}e&f\\g&h\end{bmatrix}\) \(=\) \(\begin{bmatrix} a \times e + b \times g & a \times f + b \times h \\ c \times e + d \times g & c \times f + d \times h \end{bmatrix}\)ui

這就貌似很是簡單了。spa

看一個例題,斐波那契數列不一樣於通常的斐波那契數列,\(n\)\(long\) \(long\)以內,因此\(O(n)\)絕對會超時,這是須要矩陣快速冪,複雜度是\(O(3^3logn)\) \(=\) \(O(9log n)\),忽略常數,那麼複雜度就是\(O(logn)\)code

咱們定義一個矩陣等式,而後去求問號矩陣get

\(\begin{bmatrix} f[i+1] & f[i]\end{bmatrix}\) \(\times\) \(\begin{bmatrix}? \end{bmatrix}\) \(=\) \(\begin{bmatrix}f[i+2] & f[i+1] \end{bmatrix}\)string

\(f[i+2] = f[i] + f[i+1]\)io

構造的\(?\)號矩陣就是模板

\(\begin{bmatrix} 1 & 1 \\ 1 & 0\end{bmatrix}\)class

帶回檢驗

\(\begin{bmatrix}f[i+1] & f[i] \\ 0 & 0\end{bmatrix} \times \begin{bmatrix} 1 & 1 \\ 1 & 0\end{bmatrix} = \begin{bmatrix}f[i+1] \times 1 + f[i] \times 1 & f[i+1] \times 1 + f[i] \times 0 \\ 0 \times 1 + 0 \times 1 & 0 \times 1 + 0 \times 0 \end{bmatrix}\)

上式化簡爲

\(\begin{bmatrix}f[i+2] & f[i+1] \\ 0 & 0\end{bmatrix} = \begin{bmatrix}f[i+2] & f[i+1]\end{bmatrix}\)

因此成立.

矩陣快速冪模板代碼就是

struct Mat {
    int a[3][3];
    Mat() {memset(a,0,sizeof a);}
    inline void build() {
        memset(a,0,sizeof a);
        for(re int i = 1 ; i <= 2 ; ++ i) a[i][i]=1;
    }
};
Mat operator*(Mat &a,Mat &b)
{
    Mat c;
    for(re int k = 1 ; k <= 2 ; ++ k)
        for(re int i = 1 ; i <= 2 ; ++ i)
            for(re int j = 1 ; j <= 2 ; ++ j)
                c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j]%mod)%mod;
    return c; 
}
Mat quick_Mat(int x)
{
    Mat ans;ans.build();
    while(x) {
        if((x&1)==1) ans = ans * a;
        a = a * a;
        x >>= 1;
    }
    return ans;
}

例題代碼是

#include <cstdio>
#include <iostream>
#include <cmath>
#include <cstring>
#include <queue>
#include <stack>
#define re register
#define Max 200000012
#define int long long
int n;
const int mod=1000000007;
struct Mat {
    int a[3][3];
    Mat() {memset(a,0,sizeof a);}
    inline void build() {
        memset(a,0,sizeof a);
        for(re int i = 1 ; i <= 2 ; ++ i) a[i][i]=1;
    }
};
Mat operator*(Mat &a,Mat &b)
{
    Mat c;
    for(re int k = 1 ; k <= 2 ; ++ k)
        for(re int i = 1 ; i <= 2 ; ++ i)
            for(re int j = 1 ; j <= 2 ; ++ j)
                c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j]%mod)%mod;
    return c; 
}
Mat a;
Mat quick_Mat(int x)
{
    Mat ans;ans.build();
    while(x) {
        if((x&1)==1) ans = ans * a;
        a = a * a;
        x >>= 1;
    }
    return ans;
}
signed main()
{
    scanf("%lld",&n);
    a.a[1][1]=1;a.a[1][2]=1;
    a.a[2][1]=1;Mat b;
    b.a[1][1]=1;b.a[2][1]=1;
    if(n>=1 && n<=2) {
        printf("1");return 0;
    }
    Mat ans=quick_Mat(n-2);
    ans=ans*b;
    printf("%lld",ans.a[1][1]);
    return 0;
}
相關文章
相關標籤/搜索