[SDOI 2015] 约数个数和

题目大意

d(x)d(x)xx 的约数个数,求

i=1nj=1md(ij)\sum_{i = 1}^{n} \sum_{j = 1}^{m} d(ij)

多组询问。

1n,  m,  T500,0001 \leqslant n, \; m, \; T \leqslant 500,000

题目链接

【SDOI 2015】约数个数和 - LibreOJ 2185

题解

首先有

d(nm)=injm[gcd(i,  j)=1]d(nm) = \sum_{i | n} \sum_{j | m} [gcd(i, \; j) = 1]

那么就可以开始推式子了(默认 nn 为更小的):

\begin{align} &\sum_{i = 1}^{n} \sum_{j = 1}^{m} d(ij) \\ = &\sum_{i = 1}^{n} \sum_{j = 1}^{m} \sum_{a | i} \sum_{b | j} [gcd(a, \; b) = 1] \\ = &\sum_{i = 1}^{n} \sum_{j = 1}^{m} \sum_{a | i} \sum_{b | j} \sum_{d | a, \; d | b} \mu(d) \qquad (\mu \times 1 = e) \\ = &\sum_{i = 1}^{n} \sum_{j = 1}^{m} \sum_{d | i, \; d | j} \mu(d) d(\lfloor \frac{i}{d} \rfloor) d(\lfloor \frac{j}{d} \rfloor) \\ = &\sum_{d = 1}^{n} \mu(d) \sum_{x | i} d(\lfloor \frac{i}{d} \rfloor) \sum_{d | j} d(\lfloor \frac{j}{d} \rfloor) \\ = &\sum_{d = 1}^{n} \mu(d) \sum_{i = 1}^{\lfloor \frac{n}{d} \rfloor} d(i) \sum_{j = 1}^{\lfloor \frac{m}{d} \rfloor} d(j) \end{align}

可以用线性筛在 O(n)O(n) 求出 d(x)d(x) 函数、莫比乌斯函数及其前缀和,每次询问 O(n)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
#include <cstdio>
#include <algorithm>
const int MAXN = 50005;
long long mu[MAXN], d[MAXN];
int prime[MAXN], primeCnt;
bool notPrime[MAXN];
void linearShaker() {
static int minPrimeCnt[MAXN];
notPrime[0] = notPrime[1] = true;
mu[1] = d[1] = 1;
for (int i = 2; i < MAXN; i++) {
if (!notPrime[i]) {
prime[++primeCnt] = i;
mu[i] = -1;
d[i] = 2;
minPrimeCnt[i] = 1;
}
for (int j = 1; j <= primeCnt && i * prime[j] < MAXN; j++) {
notPrime[i * prime[j]] = true;
if (i % prime[j] == 0) {
mu[i * prime[j]] = 0;
minPrimeCnt[i * prime[j]] = minPrimeCnt[i] + 1;
d[i * prime[j]] = d[i] / (minPrimeCnt[i] + 1) * (minPrimeCnt[i] + 2);
break;
} else {
mu[i * prime[j]] = -mu[i];
d[i * prime[j]] = d[i] * 2;
minPrimeCnt[i * prime[j]] = 1;
}
}
}
for (int i = 2; i < MAXN; i++) {
mu[i] += mu[i - 1];
d[i] += d[i - 1];
}
}
long long calc(int n, int m) {
if (n > m) std::swap(n, m);
long long res = 0;
for (int i = 1, last; i <= n; i = last + 1) {
last = std::min(n / (n / i), m / (m / i));
res += (mu[last] - mu[i - 1]) * d[n / i] * d[m / i];
}
return res;
}
int main() {
linearShaker();
int T;
scanf("%d", &T);
while (T--) {
int n, m;
scanf("%d %d", &n, &m);
printf("%lld\n", calc(n, m));
}
return 0;
}