HDU 5950 Recursive sequence(2016ACM/ICPC亚洲区沈阳站-重现赛（感谢东北大学）)

---

题目分析

---

这道题可以推得：$f_n = f_{n-1} + 2f_{n-2}+n^4$，因为n很大，明显可以用矩阵快速幂。但是我们会发现状态转移方程与n有关，这样我们需要利用二项式展开
$n^4 = (n-1+1)^4 = C_4^0(n-1)^4+C_4^1(n-1)^3+C_4^2(n-1)^2+C_4^3(n-1)^1+C_4^4(n-1)^0$
这样我们就可以首先构建一个初始矩阵：

$$\left\{ \begin{matrix} (n-1)^4 \\ (n-1)^3 \\ (n-1)^2 \\ (n-1)^1 \\ (n-1)^0 \\ f_{n-2} \\ f_{n-1} \\ \end{matrix} \right\} \tag{1}$$
要是初始矩阵乘以一个不变的矩阵，那么这个矩阵利用前面的二项式展开技术就可以得到了另一个用于快速幂的矩阵
$$\left\{ \begin{matrix} C_4^0& C_4^1&C_4^2&C_4^3&C_4^4&0&0 \\ 0& C_3^0&C_3^1&C_3^2&C_3^3&0&0 \\ 0& 0&C_2^0&C_2^1&C_2^2&0&0 \\ 0& 0&0&C_1^0&C_1^1&0 &0\\ 0& 0&0&0&C_0^0&0&0 \\ 0& 0&0&0&0&0&1 \\ C_4^0& C_4^1&C_4^2&C_4^3&C_4^4&2&1 \\ \end{matrix} \right\} \tag{2}$$
这样就可以矩阵快速幂了，代码写的丑，见谅！！！！

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
#define LL unsigned long long
const LL mod = 2147493647LL;

struct Matirx{
    int r, c;
    LL mat[10][10];
};

Matirx Mul(Matirx a, Matirx b){
    Matirx ret;
    ret.r = a.r;
    ret.c = b.c;
    for(int i = 0; i < ret.r; i++){
        for(int j = 0; j < ret.c; j++){
            ret.mat[i][j] = 0;
            for(int k = 0; k < a.c; k++){
                ret.mat[i][j] += a.mat[i][k]*b.mat[k][j];
                ret.mat[i][j] %= mod;
            }
        }
    }
    return ret;
}

Matirx Quick_pow(Matirx a, LL n){
    Matirx ret;
    ret.r = ret.c = 7;
    memset(ret.mat, 0, sizeof(ret.mat));
    for(int i = 0; i < 7; i++) ret.mat[i][i] = 1;
    while(n){
        if(n&1) ret = Mul(ret, a);
        a = Mul(a, a);
        n >>= 1;
    }
    return ret;
}

void build(Matirx& t){
    t.r = t.c = 7;
    t.mat[0][0] = 1, t.mat[0][1] = 4, t.mat[0][2] = 6, t.mat[0][3] = 4, t.mat[0][4] = 1, t.mat[0][5] = 0, t.mat[0][6] = 0;//1
    t.mat[1][0] = 0, t.mat[1][1] = 1, t.mat[1][2] = 3, t.mat[1][3] = 3, t.mat[1][4] = 1, t.mat[1][5] = 0, t.mat[1][6] = 0;//2
    t.mat[2][0] = 0, t.mat[2][1] = 0, t.mat[2][2] = 1, t.mat[2][3] = 2, t.mat[2][4] = 1, t.mat[2][5] = 0, t.mat[2][6] = 0;//2
    t.mat[3][0] = 0, t.mat[3][1] = 0, t.mat[3][2] = 0, t.mat[3][3] = 1, t.mat[3][4] = 1, t.mat[3][5] = 0, t.mat[3][6] = 0;//4
    t.mat[4][0] = 0, t.mat[4][1] = 0, t.mat[4][2] = 0, t.mat[4][3] = 0, t.mat[4][4] = 1, t.mat[4][5] = 0, t.mat[4][6] = 0;//5
    t.mat[5][0] = 0, t.mat[5][1] = 0, t.mat[5][2] = 0, t.mat[5][3] = 0, t.mat[5][4] = 0, t.mat[5][5] = 0, t.mat[5][6] = 1;//6
    t.mat[6][0] = 1, t.mat[6][1] = 4, t.mat[6][2] = 6, t.mat[6][3] = 4, t.mat[6][4] = 1, t.mat[6][5] = 2, t.mat[6][6] = 1;//7
}

int main(){
    int T;
    scanf("%d", &T);
    LL N, a, b;
    while(T--){
        scanf("%I64d%I64d%I64d", &N, &a, &b);
        if(N == 1) printf("%I64d\n", a%mod);
        else if(N == 2) printf("%I64d\n", b%mod);
        else{
            Matirx temp;
            build(temp);
            Matirx ret = Quick_pow(temp, N-2);
            temp.r = 7;
            temp.c = 1;
            temp.mat[0][0] = 16,temp.mat[1][0] = 8,temp.mat[2][0] = 4,temp.mat[3][0] = 2,temp.mat[4][0] = 1,temp.mat[5][0] = a%mod, temp.mat[6][0] = b%mod;
            ret = Mul(ret, temp);
            printf("%I64d\n", ret.mat[6][0]);
        }
    }
    return 0;
}