Hdu 3117 Fibonacci Numbers (线性代数_矩阵）

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=3117

题目大意：求fibonacci的第i个序列的前后各4位数，如果不多于八位数则全部输出。

解题思路：可以先求出前40个fibo数，这40个数不多于八位，全输出。如果比40大，后四位可通过矩阵快速幂求得，那前四位呢？

    以下改自http://www.cnblogs.com/zhaozhe/archive/2011/08/13/2137417.html

   Fibonacci的通项公式为：F(n)=1/sqrt(5)(((1+sqrt(5))/2)^n+((1-sqrt(5))/2)^n)

   假设F[n]可以表示成 t * 10^k（t是一个小数，k为整数），那么对于F[n]取对数log10，答案就为log10 t + K，此时很明显log10 t<1，于是我们去除整数部分，就得到了log10 t ，再用pow（10，log10 t)我们就还原回了t。将 t×1000就得到了F[n]的前四位。具体实现的时候Log10 F[n]约等于((1+sqrt(5))/2)^n/sqrt(5)，这里我们把((1-sqrt(5))/2)^n这一项忽略了，因为当N>=40时，这个数已经小的可以忽略。于是log10 F[n]就可以化简成log10 （1/sqrt(5)） + n*log10 （(1+sqrt(5))/2） = 0 - 0.5 * (log10(5.0)) + n*log10 (1+sqrt(5))/2.

   在输出的时候要注意最后四位取余完可能不足四位，要用0补上。    

测试数据：

0

40

100000000

代码：

#include <stdio.h>
#include <math.h>
#include <string.h>
#define MAX 40
#define MOD 10000
#define FIR ((1 + sqrt(5.0)) / 2)
#define SEC ((1 - sqrt(5.0)) / 2)


struct Mat{

	int mat[2][2],size;
	Mat() { 
		
		size = 2;
		mat[0][0] = 1,mat[0][1] = 1;
		mat[1][0] = 1,mat[1][1] = 0;
	}


	friend Mat operator *(Mat a,Mat b);
	friend Mat operator +(Mat a,Mat b);
	friend Mat operator ^(Mat a,int k);
}E;
Mat operator *(Mat a,Mat b){

	Mat c;
	memset(c.mat,0,sizeof(c.mat));
	for (int i = 0; i < c.size; ++i)
		for (int j = 0; j < c.size; ++j)
			for (int k = 0; k < c.size; ++k)
				if (a.mat[i][k] && b.mat[k][j])
					c.mat[i][j] = (c.mat[i][j] + a.mat[i][k] * b.mat[k][j]) % MOD;
	return c;
}
Mat operator +(Mat a,Mat b){

	Mat c;
	for (int i = 0; i < c.size; ++i)
		for (int j = 0; j < c.size; ++j)
			c.mat[i][j] = (a.mat[i][j] + b.mat[i][j]) % MOD;
	return c;
}
Mat operator ^(Mat a,int k) {

	Mat c = E;
	while (k) {

		if (k & 1) c = c * a;
		a = a * a,k >>= 1;
	}
	return c;
}


int main()
{
	int i,j,k,n;
	int last,first;
	int fibo[MAX] = {0,1};
	double fib;


	for (i = 0; i < E.size; ++i)
		E.mat[i][i] = 1,E.mat[i][1-i] = 0;
	for (i = 2; i < MAX; ++i)
		fibo[i] = fibo[i-1] + fibo[i-2];


	while (scanf("%d",&n) != EOF){

		if (n < MAX) printf("%d\n",fibo[n]);
		else {

			Mat c;
			c = c^n;
			last = c.mat[0][1] % MOD;
			//f(n)=1/sqrt(5)(((1+sqrt(5))/2)^n+((1-sqrt(5))/2)^n)
			//log10f(n) = log10 1/sqrt(5) + n*log10 (1+sqrt(5))/2
			fib = -0.5 * log10(5.0) + n * log10(FIR);
			fib = fib - (int)fib;	//log10(t)
			fib = pow(10,fib);
			while (fib < 1000) fib *= 10;
			printf("%d...",(int)fib);
			printf("%04d\n",last);
		}
	}
}

本文ZeroClock原创，但可以转载，因为我们是兄弟。