[NOI2010] 能量采集 (莫比乌斯反演 欧拉反演)

题意：

求

$$2\ast \sum _{i=1}^n\sum _{j=1}^m\gcd (i,j)-n\ast m$$

分析：

莫比乌斯反演：

那么只需要求

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

枚举每个 $\gcd (i,j)$ 的值

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

再用莫比乌斯反演，和 problem b 这题是一样的

$$\sum _{d=1}^{n}d\sum _{i=1}^{\lfloor \frac{\min (n,m)}{d}\rfloor }\mu (i)\lfloor \frac{n}{d\ast i}\rfloor \lfloor \frac{m}{d\ast i}\rfloor$$

后面可用整除分块。

代码：

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1e5 + 5;
int n, m, mobius[N], primes[N], cnt, res, sum[N];
bool st[N];
void get_mobius(int n) {
    mobius[1] = 1;
    for (int i = 2; i <= n; i ++) {
        if (!st[i]) {
            primes[cnt ++] = i;
            mobius[i] = -1;
        }
        for (int j = 0; primes[j] * i <= n; j ++) {
            int t = primes[j] * i;
            st[t] = 1;
            if (i % primes[j] == 0) {
                mobius[t] = 0;
                break;
            }
            mobius[t] = -mobius[i];
        }
    }
    for (int i = 1; i <= n; i ++) sum[i] = sum[i - 1] + mobius[i];
}
signed main() {
    get_mobius(N - 1);
    cin >> n >> m;
    for (int d = 1; d <= n; d ++) {
        int k = min(n, m) / d, ans = 0;
        for (int l = 1, r; l <= k; l = r + 1) {
            r = min(n / (n / l), m / (m / l));
            ans += (sum[r] - sum[l - 1]) * (n / (d * l)) * (m / (d * l));
        }
        res += d * ans;
    }
    cout << 2 * res - n * m << endl;
}

欧拉反演

定理：

$$n=\sum _{d\mid n}\phi (d)$$

原式为

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

将 $\gcd (i,j)$ 套用欧拉反演定理

$$\sum _{i=1}^n\sum _{j=1}^m\sum _{d\mid \gcd (i,j)}\phi (d)$$

将 $\phi (d)$ 提到前面，枚举 $d$

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

化简后面

$$\sum _{d=1}^n\phi (d)\lfloor \frac{n}{d}\rfloor \lfloor \frac{m}{d}\rfloor$$

预处理欧拉函数前缀和，用整除分块即可。

代码：

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1e5 + 5;
int n, m, primes[N], euler[N], cnt, res, sum[N];
bool st[N];
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] <= n / i; j ++) {
            int t = primes[j] * i;
            st[t] = 1;
            if (i % primes[j] == 0) {
                euler[t] = primes[j] * euler[i];
                break;
            }
            euler[t] = (primes[j] - 1) * euler[i];
        }
    }
    for (int i = 1; i <= n; i ++) sum[i] = sum[i - 1] + euler[i];
}
signed main() {
    get_eulers(N - 1);
    cin >> n >> m;
    for (int l = 1, r; l <= min(n, m); l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        res += (sum[r] - sum[l - 1]) * (n / l) * (m / l);
    }
    cout << 2 * res - n * m << endl;
}