「Luogu P2398」GCD SUM「莫比乌斯反演」

题目传送门

题意

求$\sum_{i=1}^{n}\sum_{j=1}^{n} gcd(i,j)$，$n<=10^5$

题解

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

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

$$=\sum_{d=1}^{n}d\sum_{i'=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j'=1}^{\lfloor \frac{n}{d}\rfloor} [gcd(i',j') = 1]$$

$$=\sum_{d=1}^{n}d\sum_{i'=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j'=1}^{\lfloor \frac{n}{d}\rfloor} e(gcd(i',j'))$$

$$=\sum_{d=1}^{n}d\sum_{i'=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j'=1}^{\lfloor \frac{n}{d}\rfloor} \sum_{d'|gcd(i',j')}\mu(d')$$

$$=\sum_{d=1}^{n}d\sum_{i'=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j'=1}^{\lfloor \frac{n}{d}\rfloor} \sum_{d'|i',d'|j'}\mu(d')$$

$$=\sum_{d=1}^{n} d\sum_{d'=1}^{n}\mu(d') \sum_{i'=1,d'|i'}^{\lfloor \frac{n}{d}\rfloor}\sum_{j'=1,d'|j'}^{\lfloor \frac{n}{d}\rfloor} 1$$

$$=\sum_{d=1}^{n} d\sum_{d'=1}^{n}\mu(d') \sum_{i'=1}^{\lfloor \frac{n}{dd'}\rfloor}\sum_{j'=1}^{\lfloor \frac{n}{dd'}\rfloor} 1$$

$$=\sum_{d=1}^{n} \sum_{d'=1}^{n}id(d)\mu(d') \lfloor \frac{n}{dd'}\rfloor \lfloor \frac{n}{dd'}\rfloor$$

$$=\sum_{k=1}^{n}\sum_{d=1}^{n}id(d)\mu(\lfloor \frac{k}{d}\rfloor) \lfloor \frac{n}{k}\rfloor \lfloor \frac{n}{k}\rfloor$$

$$=\sum_{k=1}^{n}(id*\mu)(k) \lfloor \frac{n}{k}\rfloor \lfloor \frac{n}{k}\rfloor$$

$$=\sum_{k=1}^{n}\phi(k) \lfloor \frac{n}{k}\rfloor \lfloor \frac{n}{k}\rfloor$$

然后可以线性筛$O(n)$预处理$\phi$，然后通过数论分块$O(\sqrt n)$计算答案.

#include <iostream>
#include <cstring>
#include <cstdio>
using namespace std;

typedef long long LL;

const int MAXN = 100100;

bool tag[MAXN];
int n, pr[MAXN], cnt, phi[MAXN];
LL phis[MAXN], ans;

void sieve(int n) { //线性筛模版 
    memset(tag, 1, sizeof tag); tag[1] = false;
    phi[1] = phis[1] = 1;
    for(int i = 2; i <= n; i ++) {
        if(tag[i]) {
            pr[++ cnt] = i;
            phi[i] = i - 1;
        }
        for(int j = 1; j <= cnt && i * pr[j] <= n; j ++) {
            tag[i * pr[j]] = false;
            if(i % pr[j] == 0) {
                phi[i * pr[j]] = phi[i] * pr[j];
                break ;
            }
            phi[i * pr[j]] = phi[i] * (pr[j] - 1);
        }
        phis[i] = phis[i - 1] + phi[i]; //phis为欧拉函数的前缀和. 
    }
}

int main() {
    scanf("%d", &n);
    sieve(n);
    for(int i = 1, j; i <= n; i = j + 1) {
        j = n / (n / i); //j为当前块的右端点，i为左端点 
        ans += (n / i) * 1ll * (n / i) * (phis[j] - phis[i - 1]);
    }
    printf("%lld\n", ans);
    return 0;
}