【luogu P3704】数字表格（数学）

数字表格

题目链接：luogu P3704

题目大意

求这个东西：
\(\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}f_{\gcd(i,j)}\)
\(f_i\) 是斐波那契额数列的第 \(i\) 位。

思路

\(\prod\limits_{i=1}^{n}\prod\limits_{j=1}^{m}f_{\gcd(i,j)}\)
\(\prod\limits_{d}f_d^{\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)==d]}\)
\(\prod\limits_{d}f_d^{\sum\limits_{i=1}^{\left\lfloor \frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{d}\right\rfloor}[\gcd(i,j)==1]}\)
\(\prod\limits_{d}f_d^{\sum\limits_{i=1}^{\left\lfloor \frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{d}\right\rfloor}\sum\limits_{d'|\gcd(i,j)}\mu(d')}\)
\(\prod\limits_{d}f_d^{\sum\limits_{d'}\mu(d')\sum\limits_{i=1}^{\left\lfloor \frac{n}{dd'}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{dd'}\right\rfloor}1}\)
\(\prod\limits_{d}f_d^{\sum\limits_{d'}\mu(d')\left\lfloor \frac{n}{dd'}\right\rfloor\left\lfloor \frac{m}{dd'}\right\rfloor}\)
这里直接搞可以搞一次询问，但是多次询问也预处理不了，考虑改改：
\(\prod\limits_{D}(\prod\limits_{d|D}f_d^{\mu(\left\lfloor \frac{D}{d}\right\rfloor)})^{\left\lfloor \frac{n}{D}\right\rfloor\left\lfloor \frac{m}{D}\right\rfloor}\)
然后括号里面的可以就预处理了。（方式是枚举 \(\left\lfloor\dfrac{D}{d}\right\rfloor\)，再枚举倍数 \(d\)）

代码

#include<cstdio>
#include<iostream>
#define ll long long
#define mo 1000000007

using namespace std;

const int N = 1e6 + 10;
int T, n, m, prime[N];
ll F[N], f[N], fv[N], mu[N];
bool np[N];

ll ksm(ll x, ll y) {
	ll re = 1;
	while (y) {
		if (y & 1) re = re * x % mo; x = x * x % mo; y >>= 1;
	}
	return re;
}

void Init() {
	f[1] = 1; fv[1] = 1; mu[1] = 1;
	for (int i = 2; i <= 1000000; i++) {
		f[i] = (f[i - 1] + f[i - 2]) % mo;
		fv[i] = ksm(f[i], mo - 2);
		if (!np[i]) mu[i] = -1, prime[++prime[0]] = i;
		for (int j = 1; j <= prime[0] && i * prime[j] <= 1000000; j++) {
			np[i * prime[j]] = 1;
			if (i % prime[j] == 0) break;
			mu[i * prime[j]] = -mu[i];
		}
	}
	F[0] = 1; for (int i = 1; i <= 1000000; i++) F[i] = 1;
	for (int i = 1; i <= 1000000; i++) {
		if (mu[i]) {
			for (int j = i; j <= 1000000; j += i)
				F[j] = F[j] * (mu[i] == 1 ? f[j / i] : fv[j / i]) % mo;
		}
		F[i] = F[i - 1] * F[i] % mo;
	}
}

ll slove(int n, int m) {
	ll re = 1;
	for (int L = 1, R; L <= n && L <= m; L = R + 1) {
		R = min(n / (n / L), m / (m / L));
		re = re * ksm(F[R] * ksm(F[L - 1], mo - 2) % mo	, 1ll * (n / L) * (m / L) % (mo - 1)) % mo;
	}
	return re;
}

int main() {
	Init();
	
	scanf("%d", &T);
	while (T--) {
		scanf("%d %d", &n, &m);
		printf("%lld\n", slove(n, m));
	}
	
	return 0;
}