[HDOJ5950]Recursive sequence（递推，二项展开，矩阵快速幂）

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=5950
题意：求解递推式f(n)=f(n-1)+2*f(n-2)+n^4。
写了个小东西，不过我的文章里式子是2*f(n-1)，内容差不多。凑合看
 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 
 4 typedef long long LL;
 5 const LL mod = 2147493647;
 6 const int maxn = 10;
 7 LL n, a, b;
 8 
 9 typedef struct Matrix {
10     LL m[maxn][maxn];
11     int r;
12     int c;
13     Matrix() {
14         r = c = 0;
15         memset(m, 0, sizeof(m));
16     }
17 } Matrix;
18 
19 Matrix mul(Matrix m1, Matrix m2) {
20     Matrix ans = Matrix();
21     ans.r = m1.r;
22     ans.c = m2.c;
23     for(int i = 1; i <= m1.r; i++) {
24         for(int j = 1; j <= m2.r; j++) {
25                for(int k = 1; k <= m2.c; k++) {
26                 if(m2.m[j][k] == 0) continue;
27                 ans.m[i][k] = ((ans.m[i][k] + m1.m[i][j] * m2.m[j][k] % mod) % mod) % mod;
28             }
29         }
30     }
31     return ans;
32 }
33 
34 Matrix quickmul(Matrix m, LL n) {
35     Matrix ans = Matrix();
36     for(int i = 1; i <= m.r; i++) {
37         ans.m[i][i]  = 1;
38     }
39     ans.r = m.r;
40     ans.c = m.c;
41     while(n) {
42         if(n & 1) {
43             ans = mul(m, ans);
44         }
45         m = mul(m, m);
46         n >>= 1;
47     }
48     return ans;
49 }
50 
51 
52 int main() {
53     // freopen("in", "r", stdin);
54     int T;
55     scanf("%d", &T);
56     while(T--) {
57         scanf("%I64d%I64d%I64d",&n,&a,&b);
58         if(n == 1) {
59             printf("%I64d
", a);
60             continue;
61         }
62         if(n == 2) {
63             printf("%I64d
", b);
64             continue;
65         }
66 
67         Matrix x; x.r = 7, x.c = 7;
68         Matrix y; y.r = 7, y.c = 1;
69         x.m[1][1]=1,x.m[1][2]=2,x.m[1][3]=1,x.m[1][4]=4,x.m[1][5]=6,x.m[1][6]=4,x.m[1][7]=1;
70         x.m[2][1]=1,x.m[2][2]=0,x.m[2][3]=0,x.m[2][4]=0,x.m[2][5]=0,x.m[2][6]=0,x.m[2][7]=0;
71         x.m[3][1]=0,x.m[3][2]=0,x.m[3][3]=1,x.m[3][4]=4,x.m[3][5]=6,x.m[3][6]=4,x.m[3][7]=1;
72         x.m[4][1]=0,x.m[4][2]=0,x.m[4][3]=0,x.m[4][4]=1,x.m[4][5]=3,x.m[4][6]=3,x.m[4][7]=1;
73         x.m[5][1]=0,x.m[5][2]=0,x.m[5][3]=0,x.m[5][4]=0,x.m[5][5]=1,x.m[5][6]=2,x.m[5][7]=1;
74         x.m[6][1]=0,x.m[6][2]=0,x.m[6][3]=0,x.m[6][4]=0,x.m[6][5]=0,x.m[6][6]=1,x.m[6][7]=1;
75         x.m[7][1]=0,x.m[7][2]=0,x.m[7][3]=0,x.m[7][4]=0,x.m[7][5]=0,x.m[7][6]=0,x.m[7][7]=1;
76 
77         y.m[1][1]=b,y.m[2][1]=a,y.m[3][1]=16,y.m[4][1]=8,y.m[5][1]=4,y.m[6][1]=2,y.m[7][1]=1;
78         Matrix p = quickmul(x,n-2);
79         Matrix ret = mul(p,y);
80         printf("%I64d
", ret.m[1][1]);
81         // cout << ret.r << " " << ret.c << endl;
82     }
83     return 0;
84 }