HDU

• 题目：
  给出n,f[1],f[2], 已知递推公式：(f[i] = 2 * f[i - 2] + f[i - 1] + i^{4});求f[n]; ((n < 2^{31}))

  通过这道题，对构造递推式的矩阵有了更深的理解。
  (f[i+1] = 2 * f[i-1] + f[i] + (i + 1)^{4}) 这个式子与上面的式子是一样了，只不过是f[i+1]。变换后我们就可以把((i+1)^{4})展开，变成 (f[i+1] = 2 * f[i-1] + f[i] + i^{4} + 4 * i^{3} + 6 * i^{2} + 4 * i + 1)。展开原因：我们要通过初始矩阵乘变换矩阵得到f[i+1]，那么最后一次乘变换矩阵前矩阵里面需要i+1前的信息，即i-1,i项的信息，不能我们求第i项,矩阵里面还要存i+1的信息，那没法推了。所以变换一下，

  将式子变成f[i+1]=第i项/第i-1项的信息，观察第3个式子，右边有七个值，所以如果要通过初始矩阵乘变换阵的到f[i+1]就要维护7个值，分别为$f[i],f[i-1],( i^{4}, i^{3}, i^{2}, i, 1 $ 将它们构成一个1x7的矩阵作为初始值阵res， 设变换阵为ex, 根据f[i]=res*ex来构造ex。ex第i列就是通过res中的各项值得到下一个res第i项的系数。以f[i+1]为例
  (f[i + 1] = f[i] + 2*f[i-1] + (i+1)^{4}) ,这里我维护第3项为((i+1)^{4}) ，ex第1列为1, 2, 1, 0, 0, 0, 0.

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<map>
#include<queue>
#include<vector>
#include<string>
#include<fstream>
using namespace std;
#define rep(i, a, n) for(int i = a; i <= n; ++ i)
#define per(i, a, n) for(int i = n; i >= a; -- i)
typedef long long ll;
const ll N = 3e12 + 105;
const double Pi = acos(- 1.0);
const int INF = 0x3f3f3f3f;
const int G = 3, Gi = 332748118;



ll mod = 2147493647;
ll qpow(ll a, ll b) { ll res = 1; while(b){ if(b) res = (res * a) % mod; a = (a * a) % mod; b >>= 1;} return res; }
ll gcd(ll a, ll b) { return b ? gcd(b, a % b) : a; }
//

struct mat{
	ll m[7][7];				//这里是矩阵大小
	mat(){memset(m,0,sizeof(m));}
};

mat mul(mat x,mat y){	//矩阵乘法
	mat res;
	for(int i=0;i<7;++i)
		for(int j=0;j<7;++j)
			for(int k=0;k<7;++k)
				res.m[i][j] = (res.m[i][j]+x.m[i][k] * y.m[k][j]) % mod;
	return res;
}

mat power(mat A, ll n){	//矩阵快速幂
	mat c = A, res;
	for(int i = 0; i < 7; ++ i) res.m[i][i]=1;	//注意res是单位阵
	while(n){
		if(n&1) res = mul(res,c);
		c = mul(c,c);
		n >>= 1;
	}
	return res;
}


ll n, A, B, T;

int main()
{
    scanf("%lld",&T);
    while(T --){
        scanf("%lld%lld%lld",&n,&A,&B);
        mat res, ex;
        res.m[0][0] = B, res.m[0][1] = A, res.m[0][2] = 81, res.m[0][3] = 27, res.m[0][4] = 9, res.m[0][5] = 3, res.m[0][6] = 1;

        ex.m[0][0] = 1, ex.m[0][1] = 1, ex.m[0][2] = 0, ex.m[0][3] = 0, ex.m[0][4] = 0, ex.m[0][5] = 0, ex.m[0][6] = 0;
        ex.m[1][0] = 2, ex.m[1][1] = 0, ex.m[1][2] = 0, ex.m[1][3] = 0, ex.m[1][4] = 0, ex.m[1][5] = 0, ex.m[1][6] = 0;
        ex.m[2][0] = 1, ex.m[2][1] = 0, ex.m[2][2] = 1, ex.m[2][3] = 0, ex.m[2][4] = 0, ex.m[2][5] = 0, ex.m[2][6] = 0;
        ex.m[3][0] = 0, ex.m[3][1] = 0, ex.m[3][2] = 4, ex.m[3][3] = 1, ex.m[3][4] = 0, ex.m[3][5] = 0, ex.m[3][6] = 0;
        ex.m[4][0] = 0, ex.m[4][1] = 0, ex.m[4][2] = 6, ex.m[4][3] = 3, ex.m[4][4] = 1, ex.m[4][5] = 0, ex.m[4][6] = 0;
        ex.m[5][0] = 0, ex.m[5][1] = 0, ex.m[5][2] = 4, ex.m[5][3] = 3, ex.m[5][4] = 2, ex.m[5][5] = 1, ex.m[5][6] = 0;
        ex.m[6][0] = 0, ex.m[6][1] = 0, ex.m[6][2] = 1, ex.m[6][3] = 1, ex.m[6][4] = 1, ex.m[6][5] = 1, ex.m[6][6] = 1;   
        
        
        if(n == 1){
            cout<<A % mod<<endl;
            continue;
        }
        else if(n == 2){
            cout<<B % mod<<endl;
            continue;
        }
        
        res = mul(res, power(ex, n - 2));
        printf("%lld
",res.m[0][0]);
    }
    return 0;
}