「Luogu P2568」GCD

Description

给定一个整数 $n\ (1 \leq n \leq 10^7)$，求满足 $\gcd(x,y)\ (1 \leq x, y \leq n)$ 为质数的数对 $(x,y)$ 的个数。

Source

[Luogu]P2568

Solution

前置知识：

常见积性函数 - 欧拉函数 - 线性筛求欧拉函数

如果 $\gcd(x,y) = 1$，那么 $\gcd(x \times p, y \times p) = p$，只需要满足 $p$ 是质数，就能符合题意。显然我们只需要找满足 $\gcd(x,y) = 1$ 的 $(x,y)$ 的对数，$\gcd(x,y) = 1$ 等价于 $x$ 与 $y$ 互质，很容易联想到 欧拉函数 。

假设 $x \leq y \leq n$，若 $\gcd(x,y) = 1$，则 $x$ 的取值有 $\varphi(y)$ 种，而 $y$ 的值可以取 $1 \sim n$，所以 $(x,y)$ 共 $\sum\limits_{i=1}^{n}\varphi(i)$ 对。在这道题中，还是假设 $x \leq y\leq n $，若 $\gcd(x \times p, y \times p) = p$，$p$ 是质数，则 $y \times p$ 一定不大于 $n$，也就是说，$y \leq \left \lfloor \frac{n}{p} \right \rfloor$ 。因为同时需要满足 $\gcd(x,y) = 1$，所以此时 $(x,y)$ 的对数为 $\sum\limits_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor}\varphi(i)$，而 $p$ 可以是 $n$ 以内的任何一个质数，因此我们需要枚举每一个 $p$，最后答案为 $\sum\limits_{p \in prime}^{n}\sum\limits_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor}\varphi(i)$ 。注意这是在 $x \leq y$ 条件下的答案，$\gcd(x,y)$ 与 $\gcd(y,x)$ 算不同的答案，计入答案时要乘 $2$，其中 $\gcd(1 \times p,1\times p)$ 被算了 $2$ 遍，所以还要减 $1$，即

$$\sum\limits_{p \in prime}^{n}\left [\left( 2\sum\limits_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor}\varphi(i)\right ) -1\right ]$$

因为要多次求 $\sum\limits_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor}\varphi(i)$，我们可以预处理出 $\varphi(1) \sim \varphi(n)$ 的 前缀和，整道题的时间复杂度为 $O(n)$ 。

Code

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;

template <class T>
inline void read(T &x) {
    x = 0;
    char c = getchar();
    bool f = 0;
    for (; !isdigit(c); c = getchar()) f ^= c == '-';
    for (; isdigit(c); c = getchar()) x = x * 10 + (c ^ 48);
    x = f ? -x : x;
}

template <class T>
inline void write(T x) {
    if (x < 0) {
        putchar('-');
        x = -x;
    }
    T y = 1;
	int len = 1;
    for (; y <= x / 10; y *= 10) ++len;
    for (; len; --len, x %= y, y /= 10) putchar(x / y + 48);
}

const int MAXN = 1e7, MAXP = 1e6;
int n, cnt, prime[MAXP], phi[MAXN];
LL ans, pre[MAXN + 5];
bool isNotPrime[MAXN + 5];

inline void getPhi(int n) {
    isNotPrime[1] = 1;
    phi[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!isNotPrime[i]) {
            prime[++cnt] = i;
            phi[i] = i - 1;
        }
        for (int j = 1; j <= cnt && i * prime[j] <= n; ++j) {
            isNotPrime[i * prime[j]] = 1;
            if (i % prime[j] == 0) {
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            } else phi[i * prime[j]] = phi[i] * (prime[j] - 1);
        }
    }
}

int main() {
    read(n);
    getPhi(n);//线性筛求 n 以内的质数 和 欧拉函数 
    for (int i = 1; i <= n; ++i) pre[i] = pre[i - 1] + phi[i];//预处理欧拉函数前缀和 
    for (int i = 1; i <= cnt; ++i)//枚举 gcd(x * prime[i], y * prime[i]) = prime[i] 
        ans += (LL) pre[n / prime[i]] * 2 - 1;//gcd(1 * prime[i], 1 * prime[i]) 算了 2 遍，要减 1 
    write(ans);
    putchar('\n');
    return 0;
}