[洛谷 P5221] Product (莫比乌斯反演)

题意

求

$$\prod _{i=1}^n\prod _{j=1}^n\frac{\text{lcm}(i,j)}{\gcd (i,j)}$$

对 $104857601$ 取模

$1\le n\le 10^6$

分析：

将 $\text{lcm}(i,j)=\frac{i\cdot j}{\gcd (i,j)}$ 替换原式得

$$\prod _{i=1}^n\prod _{j=1}^n\frac{i\cdot j}{\gcd (i,j{)}^2}$$

对 $\gcd (i,j{)}^2$ 求逆元

$$\prod _{i=1}^n\prod _{j=1}^{n}i\cdot j\cdot \gcd (i,j{)}^{-2}$$

前后分开来看

$$\prod _{i=1}^n\prod _{j=1}^{n}i\cdot j\prod _{i=1}^n\prod _{j=1}^n\gcd (i,j{)}^{-2}$$

前面 $\prod _{i=1}^n\prod _{j=1}^{n}i\cdot j=\prod _{i=1}^n{i}^n\prod _{j=1}^{n}j=\prod _{i=1}^n{i}^n\cdot n!=(n!{)}^n\prod _{i=1}^n{i}^n=(n!{)}^{2n}$

后半部分 先看原来的值再求逆元

$$\prod _{i=1}^n\prod _{j=1}^n\gcd (i,j{)}^2$$

先把平方拿出去看里面

$$\prod _{i=1}^n\prod _{j=1}^n\gcd (i,j)$$

枚举 $\gcd (i,j)$

$$\prod _{d=1}^n{d}^{\sum _{i=1}^n\sum _{j=1}^n[\gcd (i,j)=d]}$$

指数部分 因为 $\sum _{i=1}^n\sum _{i=1}^n[\gcd (i,j)=1]=2\sum _{i=1}^n\phi (i)-1$ 所以

$$\prod _{d=1}^n{d}^{2\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\phi (i)-1}$$

最后答案为

$$(n!{)}^{2n}\cdot (\prod _{d=1}^n{d}^{2\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\phi (i)-1}{)}^{-2}$$

代码(毒瘤出题人卡空间)：

#include <cstdio>
using namespace std;
const int N = 1e6 + 5, mod = 104857601;
int n, fact = 1, euler[N], primes[N], cnt, res = 1, inv = 1;
bool st[N];
int qmi(int a, int b) {
    int res = 1;
    while (b) {
        if (b & 1) res = 1ll * res * a % mod;
        a = 1ll * a * a % mod;
        b >>= 1;
    }
    return res;
}
void get_eulers(int n) {
    euler[1] = 1;
    for (int i = 2; i <= n; i ++) {
        if (!st[i]) {
            primes[cnt ++] = i;
            euler[i] = i - 1;
        }
        for (int j = 0; primes[j] * i <= n; j ++) {
            int t = primes[j] * i;
            st[t] = 1;
            if (i % primes[j] == 0) {
                euler[t] = euler[i] * primes[j];
                break;
            }
            euler[t] = euler[i] * (primes[j] - 1);
        }
    }
    for (int i = 1; i <= n; i ++) euler[i] = (euler[i - 1] + 2 * euler[i]) % (mod - 1);
}
signed main() {
    scanf("%d", &n);
    get_eulers(n);
    for (int i = 1; i <= n; i ++) fact = 1ll * fact * i % mod;
    res = qmi(fact, 2 * n);
    for (int i = 1; i <= n; i ++) {
        inv = (1ll * inv * qmi(i, euler[n / i] - 1)) % mod;
    }
    printf("%d\n", 1ll * res * qmi(1ll * inv * inv % mod, mod - 2) % mod);
}