「洛谷 P2257」YY 的 GCD - 莫比乌斯反演

发表于 分类于 ， 阅读次数： Disqus：

题意描述

洛谷链接

给定 \(N, M\)，求 \(1 \leq x \leq N\)，\(1 \leq y \leq M\) 且 \(\gcd(x, y)\) 为质数的 \((x, y)\) 有多少对。

数据范围：\(T = 10^4\)，\(N, M \leq 10^7\)。

解题思路

本题要求求下列式子的值（\(\text{Primes}\) 为质数集合）：

\[ \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) \in \text{Primes}]} } \]

看到单元函数我们可以尝试用莫比乌斯反演化简。

我们设函数 \(A(n, m) = \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) \in \text{Primes}]} }\)，进行一些化简即可得 \(A(n, m) = \sum_{ k \in \text{Primes} }{ \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) = k]} } }\)，显然 \(A(n, m)\) 即为答案。

对于 \(A(n, m)\) 这个函数的整体不好化简，我们可以先把 \(A(n, m)\) 拆出来。我们设下列函数：

\[ f(k) = \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) = k]} } \]

显然 \(A(n, m) = \sum_{k \in Primes}{f(k)}\)，故问题转化成了化简 \(f(k)\)。

我们可以尝试用莫比乌斯反演化简 \(f(k)\)。设下列函数：

\[ \begin{align*} F(k) &= \sum_{k | d}{f(d)} \\ &= \sum_{k | d}{ \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) = d]} } } \\ &= \sum_{t = 1}^{\min\{ \lfloor \frac{n}{k} \rfloor , \lfloor \frac{m}{k} \rfloor \} }{ \sum_{x = 1}^n{ \sum_{y = 1}^m{[\gcd(x, y) = tk]} } } \\ &= \sum_{x = 1}^n{ \sum_{y = 1}^m{[k | \gcd(x, y)]} } \\ &= \sum_{x = 1}^n{ \sum_{y = 1}^m{[k | x][k | y]} } \\ &= \lfloor \frac{n}{k} \rfloor \lfloor \frac{m}{k} \rfloor \end{align*} \\ \]

由于 \(F(k) = \sum_{k | d}{f(d)}\)，我们可以进行莫比乌斯反演：

\[ \begin{align*} f(k) &= \sum_{k | d}{ \mu(\frac{d}{k}) F(d) } \\ &= \sum_{k | d}{ \mu(\frac{d}{k}) \lfloor \frac{n}{d} \rfloor \lfloor \frac{m}{d} \rfloor } \\ &= \sum_{t = 1}^{ \min\{ \lfloor \frac{n}{k} \rfloor , \lfloor \frac{m}{k} \rfloor \} }{ \mu(t) \lfloor \frac{n}{tk} \rfloor \lfloor \frac{m}{tk} \rfloor } \end{align*} \]

于是我们可以得出：

\[ A(n, m) = \sum_{k \in Primes}{ \sum_{t = 1}^{ \min\{ \lfloor \frac{n}{k} \rfloor , \lfloor \frac{m}{k} \rfloor \} }{ \mu(t) \lfloor \frac{n}{tk} \rfloor \lfloor \frac{m}{tk} \rfloor } } \]

但现在直接求 \(A(n, m)\) 显然会 TLE，我们需要进行进一步优化。

我们发现在 \(A(n, m)\) 中大量出现 \(tk\)，我们可以尝试把 \(tk\) 提出来：

\[ \begin{align*} A(n, m) &= \sum_{k \in Primes}{ \sum_{t = 1}^{ \min\{ \lfloor \frac{n}{k} \rfloor , \lfloor \frac{m}{k} \rfloor \} }{ \mu(t) \lfloor \frac{n}{tk} \rfloor \lfloor \frac{m}{tk} \rfloor } } \\ &= \sum_{x = 1}^{ \min\{n, m\} }{ \sum_{ k | x, k \in \text{Primes} }{ \mu(\frac{x}{k}) \lfloor \frac{n}{x} \rfloor \lfloor \frac{m}{x} \rfloor } } \\ &= \sum_{x = 1}^{ \min\{n, m\} }{ \lfloor \frac{n}{x} \rfloor \lfloor \frac{m}{x} \rfloor \sum_{ k | x, k \in \text{Primes} }{ \mu(\frac{x}{k}) } } \end{align*} \]

对于 \(\sum_{ k | x, k \in \text{Primes} }{ \mu(\frac{x}{k}) }\)，我们可以在预处理 \(\mu(n)\) 过后预处理该式，然后对该式作前缀和。最后再用数论分块求 \(A(n, m)\) 即可。每次询问的时间复杂度为 \(O(\sqrt{n})\)。

代码演示

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
#include <cstdio>
#include <iostream>

const int MAXN = 1e7;

bool isNotPrime[MAXN + 1];
int mu[MAXN + 1], primes[MAXN + 1], cnt;
long long sumMu[MAXN + 1];

inline void getPrimes() {
    isNotPrime[0] = isNotPrime[1] = true;
    mu[1] = 1;
    for (int i = 2; i <= MAXN; i++) {
        if (!isNotPrime[i]) {
            primes[++cnt] = i;
            mu[i] = -1;
        }

        for (int j = 1; j <= cnt; j++) {
            int t = i * primes[j];
            if (t > MAXN) break;

            isNotPrime[t] = true;

            if (i % primes[j] == 0) {
                mu[t] = 0;
                break;
            } else mu[t] = -mu[i];
        }
    }

    for (int i = 1; i <= cnt; i++) {
        for (int j = 1; primes[i] * j <= MAXN; j++) {
            sumMu[primes[i] * j] += mu[j];
        }
    }
    for (int i = 2; i <= MAXN; i++) sumMu[i] += sumMu[i - 1];
}

inline long long f(int n, int m) {
    long long ans = 0;
    for (int l = 1, r = 0; l <= std::min(n, m); l = r + 1) {
        r = std::min(n / (n / l), m / (m / l));
        ans += (sumMu[r] - sumMu[l - 1]) * (n / l) * (m / l);
    }
    return ans;
}

int main() {
    getPrimes();

    int t;

    scanf("%d", &t);

    while (t--) {
        int n, m;

        scanf("%d %d", &n, &m);

        printf("%lld\n", f(n, m));
    }

    return 0;
}