Luogu5221 Product

题面

题解

首先，这种式子肯定是莫比乌斯反演之类的套路

\[\begin{aligned} & \prod_{i=1}^n\prod_{j=1}^n \frac{\mathrm{lcm}(i,j)}{\gcd(i,j)} \\ =& \prod_{i=1}^n\prod_{j=1}^n \frac{ij}{\gcd^2(i,j)} \\ \end{aligned} \]

其中

\[\begin{aligned} &\prod_{i=1}^n\prod_{j=1}^n \gcd(i,j)\\ =&\prod_{d=1}^n d^{\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d]} \end{aligned} \]

然后暴力反演或者求\(\varphi\)即可

---

这种方法太\(\texttt{naive}\)了，换一种

考虑筛出\(f(i) = \prod_{j=1}^n \gcd(i,j)\)

当\(i \in \mathrm{prime}\)时，\(f(i) = i^{n / i}\)

否则找到一个质因子\(p\)

设\(i = i' \times p^a\)，那么因为\(p\)能够加到\(\prod \gcd(i,j)\)的总贡献是\(p^{n / p^a}\)，所以

\[f(i) = f(i') \times p^{n / p^a} \]

然后暴力线性筛就可以了

时间复杂度大约是\(\mathrm{O}(n\log_2 n)\)

代码

这里面的\(ans\)就是\(f\)

#include<cstdio>
#include<cstring>
#include<cctype>
#include<algorithm>
#define RG register
#define file(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
#define clear(x, y) memset(x, y, sizeof(x))

inline int read()
{
	int data = 0, w = 1; char ch = getchar();
	while(ch != '-' && (!isdigit(ch))) ch = getchar();
	if(ch == '-') w = -1, ch = getchar();
	while(isdigit(ch)) data = data * 10 + (ch ^ 48), ch = getchar();
	return data * w;
}

const int maxn(1000010), Mod(104857601);
int prime[maxn >> 2], ans[maxn], n, mul = 1, s = 1, cnt;
bool not_prime[maxn];

int fastpow(int x, int y)
{
	int ans = 1;
	while(y)
	{
		if(y & 1) ans = 1ll * ans * x % Mod;
		x = 1ll * x * x % Mod, y >>= 1;
	}
	return ans;
}

void init()
{
	ans[1] = 1;
	for(RG int i = 2; i <= n; i++)
	{
		if(!not_prime[i]) ans[prime[++cnt] = i] = fastpow(i, n / i);
		for(RG int j = 1; j <= cnt && 1ll * i * prime[j] <= n; j++)
		{
			int k = i * prime[j], k_1 = k, p = prime[j], pow = 1;
			not_prime[k] = 1;
			while(!(k % p)) k /= p, pow *= p;
			ans[k_1] = ans[i] * fastpow(p, n / pow);
			if(!(i % prime[j])) break;
		}
	}
}

int main()
{
	n = read(); init();
	for(RG int i = 1; i <= n; i++) mul = 1ll * mul * i % Mod;
	mul = fastpow(mul, n << 1);
	for(RG int i = 2; i <= n; i++)
		s = 1ll * s * ans[i] % Mod;
	s = fastpow(s, Mod - 3);
	printf("%lld\n", 1ll * mul * s % Mod);
	return 0;
}