@bzoj - 5104@ Fib數列

目錄

• @description@
• @solution@
• @accepted code@
• @details@

---

@description@

Fib數列爲1,1,2,3,5,8...
求在Mod10^9+9的意義下，數字N在Fib數列中出如今哪一個位置
無解輸出-1

原題傳送門。

@solution@

一個熟練的 OIer 選手應該能迅速發現 5 在模 10^9 + 9 意義下有二次剩餘。考慮斐波那契通項公式：

\[f_i = \frac{1}{\sqrt{5}}((\frac{1 + \sqrt{5}}{2})^n - (\frac{1 - \sqrt{5}}{2})^n) = N \]

不妨記 \(a = \frac{1 + \sqrt{5}}{2}, b = \frac{1 - \sqrt{5}}{2}\)，則咱們要解方程 \(a^n - b^n = \sqrt{5}N = K\)。

這個方程可解嗎？注意 a, b 存在關係：\(a + b = 1, ab = -1\)。前一個顯然很差用，咱們用後一個：

\[a^n - (\frac{-1}{a})^n = K\\ a^{2n} - Ka^{n} - (-1)^n = 0\\ a^n = \frac{K \pm \sqrt{K^2 + 4\times(-1)^n}}{2} \]

討論一下 n 的奇偶能夠解出 \(a^n\)，而後 BSGS 便可。平方根也可 BSGS + 原根來解。

@accepted code@

#include <cstdio>
#include <vector>
#include <cassert>
#include <iostream>
#include <algorithm>
using namespace std;

const int MOD = int(1E9) + 9;
const int INV2 = (MOD + 1) / 2;
const int INV5 = (MOD + 1) / 5;
const int SQ5 = 383008016;
const int A = 1LL*(MOD + 1 + SQ5)*INV2%MOD;
const int B = 1LL*(MOD + 1 - SQ5)*INV2%MOD;
const int C = 1LL*SQ5*INV5%MOD;
const int G = 13;
const int BLOCK = 32000;
const int HASHSIZE = 1000037;
const int MA = 133086171;
const int INVMA = 74832817;
const int GCD = 3;

inline int add(int x, int y) {return (x + y >= MOD ? x + y - MOD : x + y);}
inline int sub(int x, int y) {return (x - y < 0 ? x - y + MOD : x - y);}
inline int mul(int x, int y) {return 1LL * x * y % MOD;}

int pow_mod(int b, int p) {
	int ret = 1;
	for(int i=p;i;i>>=1,b=mul(b,b))
		if( i & 1 ) ret = mul(ret, b);
	return ret;
}
int fib(int n) {
	int x = sub(pow_mod(A, n), pow_mod(B, n));
	return mul(C, x);
}

vector<pair<int, int> >h[HASHSIZE];
int hash_search(int x) {
	int y = x % HASHSIZE;
	for(int i=0;i<h[y].size();i++)
		if( h[y][i].first == x )
			return h[y][i].second;
	return -1;
}
void hash_insert(int x, int k) {
	int y = x % HASHSIZE;
	for(int i=0;i<h[y].size();i++)
		if( h[y][i].first == x )
			return ;
	h[y].push_back(make_pair(x, k));
}
void init() {
	int p = pow_mod(G, BLOCK);
	for(int i=1,q=p;i<=BLOCK;i++,q=mul(q,p))
		hash_insert(q, i*BLOCK);
}
int bsgs(int x) {
	if( x == 0 ) return -1;
	int ret = MOD;
	for(int i=1,q=G;i<=BLOCK;i++,q=mul(q,G)) {
		int p = hash_search(mul(x, q));
		if( p != -1 ) ret = min(ret, p - i);
	}
	assert(pow_mod(G, ret) == x);
	return ret;
}
int msqrt(int x) {
	if( x == 0 ) return 0;
	int p = bsgs(x); 
	return p % 2 ? -1 : pow_mod(G, p / 2);
}

int ans = -1;
void update(int x, int r) {
	x = 1LL*(x / GCD)*INVMA%((MOD - 1) / GCD);
	if( x % 2 == r )
		if( ans == -1 || ans > x )
			ans = x;
}
int main() {
	int N, K; scanf("%d", &N), K = mul(N, SQ5), init();
	
	int P = msqrt(add(mul(K, K), 4)); // even
	if( P != -1 ) {
		int k = bsgs(mul(add(K, P), INV2));
		if( k != -1 && k % GCD == 0 ) update(k, 0);
		k = bsgs(mul(sub(K, P), INV2));
		if( k != -1 && k % GCD == 0 ) update(k, 0);
	}
	
	P = msqrt(sub(mul(K, K), 4)); // odd
	if( P != -1 ) {
		int k = bsgs(mul(add(K, P), INV2));
		if( k != -1 && k % GCD == 0 ) update(k, 1);
		k = bsgs(mul(sub(K, P), INV2));
		if( k != -1 && k % GCD == 0 ) update(k, 1);
	}
	
	if( ans != -1 ) assert(fib(ans) == N);
	printf("%d\n", ans);
}

@details@

由於幾乎都是常數，能夠預處理出來直接用。

話說這道題的主要難點是解方程那一塊吧，感受數學味兒多一點。