HDU - 3117 Fibonacci Numbers

HDU - 3117 Fibonacci Numbers

Description

The Fibonacci sequence is the sequence of numbers such that every element is equal to the sum of the two previous elements, except for the first two elements f0 and f1 which are respectively zero and one.

What is the numerical value of the nth Fibonacci number?

Input

For each test case, a line will contain an integer i between 0 and 10 ⁸ inclusively, for which you must compute the ith Fibonacci number fi. Fibonacci numbers get large pretty quickly, so whenever the answer has more than 8 digits, output only the first and last 4 digits of the answer, separating the two parts with an ellipsis (“...”).

There is no special way to denote the end of the of the input, simply stop when the standard input terminates (after the EOF).

Sample Input

 0
1
2
3
4
5
35
36
37
38
39
40
64
65 

Sample Output

 0
1
1
2
3
5
9227465
14930352
24157817
39088169
63245986
1023...4155
1061...7723
1716...7565 

題意：求第n個斐波那契數的前4個和後4個

思路：對於前四個我們可以采用科學計數發的方式得到，

斐波那契數的通項公式是：f(n)=1/sqrt(5)(((1+sqrt(5))/2)^n+((1-sqrt(5))/2)^n)，對於40個後((1-sqrt(5))/2)^n可以忽略不計了，

後4個我們采用矩陣快速冪的方法獲得，構造的矩陣是：

#include 
#include 
#include 
#include 
#include 
using namespace std;
typedef long long ll;
const int maxn = 10;
const int mod = 10000;

int cnt;
struct Matrix {
	ll v[maxn][maxn];
	Matrix() {}
	Matrix(int x) {
		init();
		for (int i = 0; i < maxn; i++) 
			v[i][i] = x;
	}
	void init() {
		memset(v, 0, sizeof(v));
	}
	Matrix operator *(Matrix const &b) const {
		Matrix c;
		c.init();
		for (int i = 0; i < cnt; i++)
			for (int j = 0; j < cnt; j++)
				for (int k = 0; k < cnt; k++)
					c.v[i][j] = (c.v[i][j] + (ll)v[i][k]*b.v[k][j]) % mod;
		return c;
	}
	Matrix operator ^(int b) {
		Matrix a = *this, res(1);
		while (b) {
			if (b & 1)
				res = res * a;
			a = a * a;
			b >>= 1;
		}
		return res;
	}
} a, b, tmp;

int main() {
	int f[40], n;
	f[0] = 0, f[1] = 1;
	for (int i = 2; i < 40; i++) 
		f[i] = f[i-1] + f[i-2];
	while (scanf("%d", &n) != EOF) {
		if (n < 40) {
			printf("%d\n", f[n]);
			continue;
		}
		double k = log10(1.0/sqrt(5.0)) + (double)n * log10((1.0 + sqrt(5.0))/2.0);
	 	double num = k;
		num = k - (int)num;
		a.init();
		a.v[0][0] = a.v[0][1] = a.v[1][0] = 1;
		cnt = 2;
		tmp = a^n;
		printf("%d...%0.4d\n", (int)(1000.0*pow(10.0, num)), tmp.v[0][1]%mod);
	}
}