P4000 斐波那契数列-526互联

题意

求

\[fib_n \pmod{p} \]

\[n \leqslant 10^{30000000}, p < 2^{31} \]

Solution

idea by Itst

对于如此庞大的 \(n\) ，一个显而易见的思路为寻找斐波那契数列在模数 \(p\) 下的循环节。

我们可以做以下简便计算，如果说二元组 \((fib_n, fib_{n+1})\) 与 \((fib_m, fib_{m+1})\) ，满足 \(fib_n = fib_m, fib_{n+1} = fib_{m+1}\)，我们可以认为 \(\left\vert n - m \right\vert\) 为斐波那契数列数列在模数 \(p\) 下的循环节。

根据斐波那契数列的定义，显然这个简便运算的置信度是 \(100\%\) ，因为后续的数列完全可以由这两个数递推产生，结果应该是完全一模一样的。

不难考虑可以使用map，维护二元组的值，在 \(O(\log{p})\) 的时间复杂度下进行比较，考虑到 \(p\) 的范围，可以使用一个64位数来表示二元组。

根据生日悖论，随着随机数量的增加，命中率呈平方级增长，上述随机的期望时间复杂度为 \(O(\sqrt{p})\)。

现在来计算一下总体的时间复杂度，大致为 \(O(\sqrt{p} · \log^2{p} · 2^3)\)，大约为 \(2e9\) 级别，寄咯！

提供一种 \(O(\sqrt{p})\) 预处理， \(O(1)\) 查询的快速幂方法 LOJ 快速幂 2，便可以处理掉此题。

Code

点击查看代码

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <random>
#include <unordered_map>

using namespace std;

const int N = 3e7 + 50, S = 3e5 + 50;
const int mod = 998244352;
typedef long long lld;
typedef unsigned long long ulld;

inline int read() {
	register int w = 0, f = 1;
	register char c = getchar();
	while (c > '9' || c < '0') {
		if (c == '-')  f = -1;
		c = getchar();
	}
	while (c >= '0' && c <= '9') {
		w = w * 10 + c - '0';
		c = getchar();
	}
	return w * f;
}

struct Mat {
	int n, m;
	int dat[3][3];
	Mat () {
		n = m = 0;
		memset(dat, 0, sizeof (dat));
	}
	int* operator [] (register int x)  {
		return dat[x];
	}
};

int p;

unordered_map <ulld, lld> mp;

Mat operator * (register Mat a, register Mat b) {
	Mat c;
	c.n = a.n, c.m = b.m;
	for (register int i = 1; i <= c.n; ++i)
		for (register int j = 1; j <= c.m; ++j)
			for (register int k = 1; k <= a.m; ++k)
				c[i][j] = (c[i][j] + 1ll * a[i][k] * b[k][j] % p) % p;
	return c;
}

mt19937_64 rnd(time(0));
char s[N];

int Sp;
lld len, sum;

Mat A[S], B[S];
Mat Base, fib;

int main() {
	scanf("%s %d", s + 1, &p);
	Sp = (1 << 18);
	fib.n = 1, fib.m = 2;
	fib[1][1] = fib[1][2] = 1;
	A[0][1][1] = A[0][2][2] = 1;
	B[0][1][1] = B[0][2][2] = 1;
	Base[1][1] = Base[1][2] = Base[2][1] = 1;
	Base.n = Base.m = A[0].n = A[0].m = B[0].n = B[0].m = 2;
	for (register int i = 1; i <= Sp; ++i)  A[i].n = A[i].m = 2, A[i] = A[i - 1] * Base;
	for (register int i = 1; i <= Sp; ++i)  B[i].n = B[i].m = 2, B[i] = B[i - 1] * A[Sp];

	while (1) {
		register lld x = ((rnd() << 28ll) >> 28ll);	
		register Mat C = A[x & (Sp - 1)] * B[x >> 18];
		C = fib * C;
		register ulld val = ((1ull * C[1][1]) << 32) | C[1][2];
		if (mp.find(val) != mp.end()) {
			len = abs(mp[val] - x);
			break;
		}
		mp[val] = x;
	}

	for (register int i = 1; s[i]; i++)
		sum = (sum * 10 + s[i] - '0') % len;
	sum -= 2;
	printf("%d\n", (fib * (A[sum & (Sp - 1)] * B[sum >> 18]))[1][1]);
	return 0;
}