hdu 5950 Recursive sequence 矩阵快速幂

题目链接：

http://acm.hdu.edu.cn/showproblem.php?pid=5950

题意：

已知递推公式：F(n) = 2*F(n-2) + F(n-1) + n⁴和F(1) = a，F(2) = b；给定一个N，求F(N)等于多少？

思路：

我xjb凑的，7*7的矩阵，主要是n^4怎么拆，这里有个很有道理的拆开n^4：http://blog.csdn.net/spring371327/article/details/52973534

代码：

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 typedef long long ll;
 4 #define MS(a) memset(a,0,sizeof(a))
 5 #define MP make_pair
 6 #define PB push_back
 7 const int INF = 0x3f3f3f3f;
 8 const ll INFLL = 0x3f3f3f3f3f3f3f3fLL;
 9 inline ll read(){
10     ll x=0,f=1;char ch=getchar();
11     while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
12     while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
13     return x*f;
14 }
15 //////////////////////////////////////////////////////////////////////////
16 const int maxn = 1e5+10;
17 const ll mod = 2147493647;
18 
19 ll N,A,B;
20 
21 struct node{
22     ll m[8][8];
23     void clear(){
24         MS(m);
25     }
26 }a;
27 
28 node mul(node a,node b){
29     node re; re.clear();
30     for(int i=0; i<7; i++){
31         for(int k=0; k<7; k++){
32             if(a.m[i][k]){
33                 for(int j=0; j<7; j++)
34                     re.m[i][j] = (re.m[i][j]+(a.m[i][k]*b.m[k][j])%mod)%mod;
35             }
36         }
37     }
38     return re;
39 }
40 
41 node qpow(ll k){
42     node re;
43     re.clear();
44     for(int i=0; i<7; i++) re.m[i][i]=1;
45     while(k){
46         if(k&1) re=mul(re,a);
47         a = mul(a,a);
48         k >>= 1;
49     }
50     return re;
51 }
52 
53 int main(){
54     int T = read();
55     while(T--){
56         a.clear();
57         cin>>N>>A>>B;
58         a.m[0][0]=1,a.m[0][1]=2,a.m[0][2]=1,a.m[0][3]=8,a.m[0][4]=24,a.m[0][5]=32,a.m[0][6]=16;
59         a.m[1][0]=1,a.m[1][1]=0,a.m[1][2]=0,a.m[1][3]=0,a.m[1][4]=0,a.m[1][5]=0,a.m[1][6]=0;
60         a.m[2][0]=0,a.m[2][1]=0,a.m[2][2]=1,a.m[2][3]=4,a.m[2][4]=6,a.m[2][5]=4,a.m[2][6]=1;
61         a.m[3][0]=0,a.m[3][1]=0,a.m[3][2]=0,a.m[3][3]=1,a.m[3][4]=3,a.m[3][5]=3,a.m[3][6]=1;
62         a.m[4][0]=0,a.m[4][1]=0,a.m[4][2]=0,a.m[4][3]=0,a.m[4][4]=1,a.m[4][5]=2,a.m[4][6]=1;
63         a.m[5][0]=0,a.m[5][1]=0,a.m[5][2]=0,a.m[5][3]=0,a.m[5][4]=0,a.m[5][5]=1,a.m[5][6]=1;
64         a.m[6][0]=0,a.m[6][1]=0,a.m[6][2]=0,a.m[6][3]=0,a.m[6][4]=0,a.m[6][5]=0,a.m[6][6]=1;
65         ll r1=B, r2=A;
66         if(N==1) { cout<<r2<<endl; }
67         else if(N==2) { cout<<r1<<endl; }
68         else{
69             a = qpow(N-2);
70             // cout << a.m[0][0] << " " << a.m[0][1] << " "<< a.m[0][2] << " " << a.m[0][3] << " "<< a.m[0][4] << " " << a.m[0][5] << " " << a.m[0][6] << endl;
71             ll ans = ((a.m[0][0]*B)%mod+(a.m[0][1]*A)%mod+a.m[0][2]+a.m[0][3]+a.m[0][4]+a.m[0][5]+a.m[0][6])%mod;
72             cout << ans << endl;
73         }
74     }
75 
76     return 0;
77 }