HDU5950 矩陣快速冪（巧妙的遞推）

題目連接：http://acm.hdu.edu.cn/showproblem.php?pid=5950

題意：f[n] = 2*f[n-2] + f[n-1] + n^4

思路：對於遞推題而言，若是遞推n次很大，則考慮矩陣快速冪的方式推出遞推式，計算出累乘的矩陣

本題遞推式：本題的遞推式子雖然已經給出，可是因爲n^4的關係，直接是沒法使用這個f[n] = 2*f[n-2] + f[n-1] + n^4遞推完成矩陣的推導的，而是能夠先處理一下，以下：

f[n] = 2*f[n-2] + f[n-1] + n^4

f[n+1] = 2*f[n-1] + f[n] + n^4 + 4*n^3 + 6*n^2 + 4*n + 1

f[n+2] = 2*f[n] + f[n+1] + (n+1)^4 + 4*(n+1)^3 + 6*(n+1)^2 + 4*(n+1)+ 1

此時，咱們發現從n+1項開始包括n+1項，都是由7個部分組成的多項式，則咱們能夠利用n+1項和n+2項的多項式進行矩陣快速冪的遞推矩陣的推導，因爲矩陣乘法的性質，對於一個1X7的矩陣A，要求相乘另外一個矩陣B以後，仍是一個1X7的矩陣，則矩陣B的規模必須是7X7，下面是推導， 對於黃色的一行乘綠色一列，獲得橙色的一個數

 

 

完成矩陣的遞推以後，就很簡單了，用矩陣的快速冪計算便可，須要注意的是對於n>=3，咱們才須要進行矩陣相乘的運算，而初始的時候，咱們須要計算出黃色矩陣表明的部分，本題就是將n==2代入，算出初始黃色矩陣爲[a, b, 16, 8, 4, 2, 1]

代碼：

 1 #include<iostream>
 2 #include<stdio.h>
 3 using namespace std;
 4 
 5 const long long mod = 2147493647;
 6 struct mat{
 7     long long m[7][7];
 8 };
 9 
10 mat operator * (mat a, mat b){        //重載乘號，同時將數據mod10000 
11     mat ret;
12     for(int i = 0; i < 7; i++){
13         for(int j = 0; j < 7; j++){
14             long long temp = 0;
15             for(int k = 0; k < 7; k++){
16                 temp += a.m[i][k] * b.m[k][j];
17                 temp %= mod;
18             }
19             ret.m[i][j] = temp;
20         }
21     }    
22     return ret;
23 }
24 
25 mat pow_mat(int f1, int f2, mat a, int n){        //矩陣快速冪和快速冪相同（廣義快速冪的思想） 
26     mat res;
27     res.m[0][0] = f1,res.m[0][1] = f2,res.m[0][2] = 16,res.m[0][3] = 8,res.m[0][4] = 4,res.m[0][5] = 2,res.m[0][6] = 1;
28     while(n){
29         if(n&1) res = res * a;
30         a = a*a;
31         n >>= 1;
32     }
33     return res;
34 }
35 
36 int main(){
37     int t;
38     scanf("%d", &t);
39     while(t--){
40         int n, a, b;
41         scanf("%d%d%d", &n, &a, &b);
42         if(n == 1) printf("%d\n", a);
43         else if(n == 2) printf("%d\n", b);
44         else{
45             mat x;
46             x.m[0][0] = 0,x.m[0][1] = 2,x.m[0][2] = 0,x.m[0][3] = 0,x.m[0][4] = 0,x.m[0][5] = 0,x.m[0][6] = 0;
47             x.m[1][0] = 1,x.m[1][1] = 1,x.m[1][2] = 0,x.m[1][3] = 0,x.m[1][4] = 0,x.m[1][5] = 0,x.m[1][6] = 0;
48             x.m[2][0] = 0,x.m[2][1] = 1,x.m[2][2] = 1,x.m[2][3] = 0,x.m[2][4] = 0,x.m[2][5] = 0,x.m[2][6] = 0;
49             x.m[3][0] = 0,x.m[3][1] = 4,x.m[3][2] = 4,x.m[3][3] = 1,x.m[3][4] = 0,x.m[3][5] = 0,x.m[3][6] = 0;
50             x.m[4][0] = 0,x.m[4][1] = 6,x.m[4][2] = 6,x.m[4][3] = 3,x.m[4][4] = 1,x.m[4][5] = 0,x.m[4][6] = 0;
51             x.m[5][0] = 0,x.m[5][1] = 4,x.m[5][2] = 4,x.m[5][3] = 3,x.m[5][4] = 2,x.m[5][5] = 1,x.m[5][6] = 0;
52             x.m[6][0] = 0,x.m[6][1] = 1,x.m[6][2] = 1,x.m[6][3] = 1,x.m[6][4] = 1,x.m[6][5] = 1,x.m[6][6] = 1;
53             mat ans = pow_mat(a, b, x, n-2);
54             printf("%d\n", ans.m[0][1]);
55         }
56     }
57     return 0;
58 }