HDU 2030 Happy Necklace（矩阵快速幂

题目来源：http://acm.hdu.edu.cn/showproblem.php?pid=6030

题目描述：

串珠子，非环形，每个素数的连续子序列，红珠子的数量不小于蓝珠子的数量即可算一种，问给定长度为n能有多少种方案，结果对 $10^9+7$109+7取余

思路：

数据大，一看就是矩阵快速幂，推出前５项，F(2)=3,F(3)=4,F(4)=6,F(5)=9,F(6)=13....得到
F(n)=F(n-1)+F(n-3)从而得到中间矩阵A, 答案(F(n))就是 Matrix ans = $A^{n-4}\times B$An−4×B的第一项 :ans.m[0][0]

$\because \left\{\begin{array}{c}F(n-1)\\ F(n-2)\\ F(n-3)\\ \dots \end{array}\right\}\times A=\left\{\begin{array}{c}F\left(n\right)\\ F(n-1)\\ F(n-2)\\ \dots \end{array}\right\}$∵⎩⎪⎪⎨⎪⎪⎧​F(n−1)F(n−2)F(n−3)…​⎭⎪⎪⎬⎪⎪⎫​×A=⎩⎪⎪⎨⎪⎪⎧​F(n)F(n−1)F(n−2)…​⎭⎪⎪⎬⎪⎪⎫​
$又\because \left\{\begin{array}{c}F(n-1)\\ F(n-2)\\ F(n-3)\end{array}\right\}\times \left\{\begin{array}{ccc}1& 0& 1\\ 1& 0& 0\\ 0& 1& 0\end{array}\right\}=\left\{\begin{array}{c}F\left(n\right)\\ F(n-1)\\ F(n-2)\end{array}\right\}$又∵⎩⎨⎧​F(n−1)F(n−2)F(n−3)​⎭⎬⎫​×⎩⎨⎧​110​001​100​⎭⎬⎫​=⎩⎨⎧​F(n)F(n−1)F(n−2)​⎭⎬⎫​
$\therefore A=\left\{\begin{array}{ccc}1& 0& 1\\ 1& 0& 0\\ 0& 1& 0\end{array}\right\}$∴A=⎩⎨⎧​110​001​100​⎭⎬⎫​
$B=\left\{\begin{array}{c}F\left(4\right)=6\\ F\left(3\right)=4\\ F\left(2\right)=3\end{array}\right\}$B=⎩⎨⎧​F(4)=6F(3)=4F(2)=3​⎭⎬⎫​

参考代码：

#include<cstdio>
#include<cstring>
#include<iostream>
using namespace std;
#define ll long long
#define endl '\n'
#define mem(a,b) memset(a,b,sizeof(a))
const ll N=3,mod=1e9+7,INF=1e9;
struct Matrix {
    ll m[N][N];

    Matrix() {
        mem(m, 0);
    }

    void init() {
        for (int i = 0; i < N; i++)m[i][i] = 1;
    }

    Matrix operator+(const Matrix &b) const {
        Matrix c;
        for (int i = 0; i < N; i++) {
            for (int j = 0; j < N; j++) {
                c.m[i][j] = (m[i][j] + b.m[i][j]) % mod;
            }
        }
        return c;
    }

    Matrix operator*(const Matrix &b) const {
        Matrix c;
        mem(c.m,0);
        for (int i = 0; i < N; i++) {
            for (int j = 0; j < N; j++) {
                for (int k = 0; k < N; k++) {
                    c.m[i][j] = (c.m[i][j] + (m[i][k] * b.m[k][j]) % mod) % mod;
                }
            }
        }
        return c;
    }

    Matrix operator^(const ll &t) const {
        Matrix ans, a = (*this);
        ans.init();
        ll n = t;
        while (n) {
            if (n & 1) ans = ans * a;
            a = a * a;
            n >>= 1;
        }
        return ans;
    }
};
int main() {
    int t;
    cin >> t;
    while (t--) {
        ll n;
        cin>>n;
        if (n==2){
            cout<<"3"<<endl;
        }
        else if (n==3){
            cout<<"4"<<endl;
        }
        else if (n==4){
            cout<<"6"<<endl;
        }
        else {
            Matrix A;
            int at[10][10] = {{1, 0, 1},
                              {1, 0, 0},
                              {0, 1, 0}};
            for (int i = 0; i < 3; i++)
                for (int j = 0; j < 3; j++)
                    A.m[i][j] = at[i][j];
            Matrix B;
            mem(B.m,0);
            B.m[0][0] = 6;
            B.m[1][0] = 4;
            B.m[2][0] = 3;
            Matrix ans = (A ^ (n - 4)) * B;
            cout << ans.m[0][0] % mod << endl;
        }
    }
    return 0;
}

下面是一个大佬的做法：https://blog.csdn.net/survivorone/article/details/71375092

很是神奇值得学习。