「GCJ 2008 Round 1A C」numbers

Description

计算实数 ((3 + sqrt{5}) ^n) 的整数最后三位并输出，必要是补足前导 0。

Hint

Small

• (2le nle 30)

Large

• (2le nle 2 imes 10^9)

Solution

显然不可以直接算。

---

• 当 (n=1) 时， ((3 +sqrt{5})^1 = 3+sqrt{5})
• 当 (n=2) 时， ((3 +sqrt{5})^2 = 14+6sqrt{5})
• 当 (n=3) 时， ((3 +sqrt{5})^3 = 82+32sqrt{5})
• (cdots cdots)

可以发现，((3+5sqrt{5})^n) 可以表示为 (a_n + b_nsqrt{5}) 的形式。

同理， ((3-sqrt{5})^n) 也可以表示为 (a_n - b_nsqrt{5}) 的形式。

以上 (a_n,b_n) 均为整数。

那么把两者相加，得 ((3+sqrt{5})^n + (3-sqrt{5})^n = a_n + b_nsqrt{5} + a_n - b_nsqrt{5} = 2a_n) 。

又因为 (3 - sqrt{5} approx ‭0.763 < 1) ，那么 ((3-sqrt{5})^n<1) 。

由 ((3+sqrt{5})^n = 2a_n - (3-sqrt{5})^n) 不难得到 ((3+sqrt{5})^n) 的整数部分为 (2a_n - 1) 。

---

显然 ((3 + sqrt{5})^n) 可以由 ((3+sqrt{5})^{n-1}) 乘上一个 ((3+sqrt{5})) 得到。

那么： ((3+sqrt{5})^{n-1}(3+sqrt{5}) = (a_{n-1}+b{n-1}sqrt{5})(3+sqrt{5}) = (3a_{n-1} + 5b_{n-1}) + (a_{n-1}+3b_{n-1}) sqrt{5})

于是： (egin{cases} a_n = 3a_{n-1} + 5b_{n-1} \ b_n = a_{n-1}+3b_{n-1} end{cases})

这就是个递推式子，答案就是 ((2a_{n-1} -1 ) mod 10^3) 的值。

---

large 数据 (Theta(n)) 递推绝对爆炸，怎么办？

不慌，我们可以 矩阵加速！

由上面递推式可得：

[egin{bmatrix} a_n \ b_n end{bmatrix} = egin{bmatrix} 3 & 5 \ 1 & 3 end{bmatrix} imes egin{bmatrix} a_{n-1} \ b_{n-1}end{bmatrix} ]

时间复杂度直降至 (Theta (log n)) ，可以通过。

Code

/*
 * Author : _Wallace_
 * Source : https://www.cnblogs.com/-Wallace-/
 * Problem : GCJ 2008 Round 1A C numbers
 */
#include <iostream>
#include <cstdio>
#include <vector>
using namespace std;

const int mod = 1e3;
typedef vector<int> vec;
typedef vector<vec> mat;

inline mat operator * (mat a,mat b) {
	mat x(2, vec(2, 0));
	for (register int i = 0; i < 2; i++)
		for (register int j = 0; j < 2; j++)
			for (register int k = 0; k < 2; k++)
				(x[i][j] += a[i][k] * b[k][j]) %= mod;
	return x;
}

inline mat matI() {
	mat x(2, vec(2, 0));
	x[0][0] = x[1][1] = 1;
	return x;
}

mat fastpow(mat base, int k) {
	if (!k) return matI();
	mat x = fastpow(base, k >> 1);
	if (k & 1) return (x * x) * base;
	else return x * x;
}

signed main() {
	int n, T;
	cin >> T;
	for (register int i = 1; i <= T; i++) {
		cin >> n;
		mat base(2, vec(2, 0));
		base[0][0] = 3, base[0][1] = 5;
		base[1][0] = 1, base[1][1] = 3;
                /*因为 a_1 = b_1 = 1 ，所以一开始的矩阵乘不乘无所谓*/
		base = fastpow(base, n);
		printf("Case #%d: %03d
", i, (base[0][0] * 2 - 1) % mod);
	}
}