矩陣相乘，快速演算法 HDOJ 4291 A Short problem

對於g（x） % mod；很明顯，x一定有一個迴圈結；

證明：

只要存在g（x）= g（k）， g（x - 1） =  g（k - 1），那麼就會出現迴圈，而g（x）的值是在有限域（1 - mod）內變化，那麼g(x) ，g(x - 1)的排列也是有限的，而x是無限的，所以一定會出現g（x）= g（k）， g（x - 1） =  g（k - 1）；

 

那麼可以求出每層函數中的mod，依次為

mod1 = 1000000007;

mod2 = 222222224；
mod3 = 183120;
mod4 = 240;

所以i層函數中餘modi，中間值餘mod（i-1）;

 

而g（n） = 3 * g(n - 1) + g(n - 2);那麼有

g[1]                 1                  g[n]                          3          1                            g[n - 1]

             =                                                   =                                      *             

g[0]                 0                  g[n - 1]                     1           0                          g[n - 2]

 

可以得到矩陣的冪，類似與快速冪計算方法，二分計算

 

#include<stdio.h>const __int64 mod1 = 1000000007;const __int64 mod3 = 183120;const __int64 mod2 = 222222224, mod4 = 240;__int64 result[241];int copy(__int64 b[][2], __int64 a[][2]){b[0][0] = a[0][0];b[0][1] = a[0][1];b[1][0] = a[1][0];b[1][1] = a[1][1];return 0;}__int64 f(__int64 a[][2], __int64 x, __int64 b[][2], __int64 mod){__int64 temp[2][2];if(x == 1){copy(b, a);return 0;}f(a, x / 2, temp, mod);b[0][0] = (temp[0][0] * temp[0][0] + temp[0][1] * temp[1][0]) % mod;b[0][1] = (temp[0][0] * temp[0][1] + temp[0][1] * temp[1][1]) % mod;b[1][0] = (temp[1][0] * temp[0][0] + temp[1][1] * temp[1][0]) % mod;b[1][1] = (temp[1][0] * temp[0][1] + temp[1][1] * temp[1][1]) % mod;if(x % 2 == 0)return 0;copy(temp, b);b[0][0] = (temp[0][0] * a[0][0] + temp[0][1] * a[1][0]) % mod;b[0][1] = (temp[0][0] * a[0][1] + temp[0][1] * a[1][1]) % mod;b[1][0] = (temp[1][0] * a[0][0] + temp[1][1] * a[1][0]) % mod;b[1][1] = (temp[1][0] * a[0][1] + temp[1][1] * a[1][1]) % mod;return 0;}__int64 g(__int64 x, __int64 mod){__int64 a[2][2], b[2][2];if(x < 2)return x;a[0][0] = 3;a[0][1] = 1;a[1][0] = 1;a[1][1] = 0;f(a, x - 1, b, mod);return b[0][0];}int init(){int i;result[0] = 0;result[1] = 1;for(i = 2; i <= 240; i ++){result[i] = g(g(g(i % mod4, mod3) % mod3, mod2) % mod2, mod1) % mod1;}return 0;}int main(){__int64 n;init();while(scanf("%I64d", &n) != EOF){printf("%I64d\n", result[n % mod4]);}return 0;}