[SDOI2015] 约数个数和 (莫比乌斯反演)

题意：

设 $d(x)$ 为 $x$ 的约数个数，求

$$\sum _{i=1}^N\sum _{j=1}^{M}d(ij)$$

分析：

首先 $i=p_1^{{\alpha }_1}\cdots p_k^{{\alpha }_k},j=p_1^{{\beta }_1}\cdots p_k^{{\beta }_k},i\times j=p_1^{{\alpha }_1+{\beta }_1}\cdots p_k^{{\alpha }_k+{\beta }_k}$

约数个数和为

$$d(ij)=({\alpha }_1+{\beta }_1+1)\cdots ({\alpha }_k+{\beta }_k+1)={\alpha }_1{\beta }_1+{\alpha }_1{\beta }_2+\cdots +1$$

所以$$d(ij)=\sum _{x|i}\sum _{y|j}[\gcd (x,y)=1]$$

所以原式为

$$\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[\gcd (x,y)=1]$$

设 $g(n)$ 为：

$$g(n)=\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[n\mid \gcd (x,y)]$$

设 $f(n)$ 为：

$$\sum _{i=1}^N\sum _{j=1}^M\sum _{x|i}\sum _{y|j}[\gcd (x,y)=n]$$

则有

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

那么就可以莫比乌斯反演了

$$f(n)=\sum _{n\mid d}\mu (\frac{d}{n})g(d)$$

交换 $g(n)$ 的求和次序

由于 $x,y$ 分别是 $i,j$ 的约数，所以可以先枚举 $x,y$

那么后面的 $[n\mid \gcd (x,y)]$ 与 $i,j$ 无关不需要考虑，所以只需要计算 $i$ 中有多少 $x$ 的倍数， $j$ 中有多少 $y$ 的倍数，所以是 $\lfloor \frac{N}{x}\rfloor \lfloor \frac{M}{y}\rfloor$

所以式子就变为了

$$\sum _{x=1}^N\sum _{y=1}^M\lfloor \frac{N}{x}\rfloor \lfloor \frac{M}{y}\rfloor [n\mid \gcd (x,y)]$$

$[n\mid \gcd (x,y)]$ 只需要关心 $n$ 的倍数即可，那么就是

$$x^{\prime }=\frac{x}{n},y^{\prime }=\frac{y}{n}$$

替换得

$$\sum _{x^{\prime }=1}^{\frac{N}{n}}\sum _{y^{\prime }=1}^{\frac{M}{n}}\lfloor \frac{N}{n{x}^{\prime }}\rfloor \lfloor \frac{M}{n{y}^{\prime }}\rfloor$$

设

$$N^{\prime }=\frac{N}{n},M^{\prime }=\frac{M}{n}$$

则

$$\sum _{x^{\prime }=1}^{N^{\prime }}\sum _{y^{\prime }=1}^{M^{\prime }}\lfloor \frac{N^{\prime }}{x^{\prime }}\rfloor \lfloor \frac{M^{\prime }}{y^{\prime }}\rfloor$$

假设有二重积分

$${\iint }_{D}f(x,y)\text{d}x\text{d}y$$

当区域 $D$ 为矩形区域时，可以转为二次积分

$${\int }_{D_1}f(x,y)\text{d}x{\int }_{D_2}f(x,y)\text{d}y$$

那么原式可以变为

$$\sum _{x^{\prime }=1}^{N^{\prime }}\lfloor \frac{N^{\prime }}{x^{\prime }}\rfloor \sum _{y^{\prime }=1}^{M^{\prime }}\lfloor \frac{M^{\prime }}{y^{\prime }}\rfloor$$

设

$$h(x)=\sum _{i=1}^x\lfloor \frac{x}{i}\rfloor$$

那么

$$f(n)=\sum _{n\mid d}\mu (\frac{d}{n})g(d)$$

答案为

$$f(1)=\sum _{d=1}^N\mu (d)g(d)$$

带入 $g(d)$

$$\sum _{i=1}^{\min (N,M)}\mu (i)h(\lfloor \frac{N}{i}\rfloor )h(\lfloor \frac{M}{i}\rfloor )$$

可以用整除分块计算，$h(x)$ 也同样可以用整除分块预处理。

代码：

#include <stdio.h>
#include <algorithm>
#define int long long
using namespace std;
const int N = 5e4 + 5;
int T, n, m, primes[N], mobius[N], cnt, sum[N], h[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];
}
void get_h(int n) {
    for (int i = 1; i <= n; i ++) {
        for (int l = 1, r; l <= i; l = r + 1) {
            r = min(i, i / (i / l));
            h[i] += (r - l + 1) * (i / l);
        }
    }
}
signed main() {
    get_mobius(N - 1), get_h(N - 1);
    scanf("%lld", &T);
    while (T --) {
        scanf("%lld%lld", &n, &m);
        int res = 0, k = min(n, m);
        for (int l = 1, r; l <= k; l = r + 1) {
            r = min(k, min(n / (n / l), m / (m / l)));
            res += (sum[r] - sum[l - 1]) * h[n / l] * h[m / l];
        }
        printf("%lld\n", res);
    }
}