求
$$\sum_{i = 1} ^ n\sum_{j = 1} ^ nd(ij)$$
链接
BZOJ 4176
题解
首先有一个结论
$$d(nm) = \sum_{i | n}\sum_{j | m}[\gcd(i, j) = 1]$$
证明：设 $nm = p_1^{x_1}\cdot p_2^{x_2}\cdot p_3^{x_3}\cdots p_k^{x_k}$，$n = p_1^{y_1}\cdot p_2^{y_2}\cdot p_3^{y_3}\cdots p_k^{y_k}$，则
$m = p_1^{x_1 - y_1}\cdot p_2^{x_2 - y_2}\cdot p_3^{x_3 - y_3}\cdots p_k^{x_k - y_k}$。
设 $i = p_1^{a_1}\cdot p_2^{a_2}\cdot p_3^{a_3}\cdots p_k^{a_k}$，$j = p_1^{b_1}\cdot p_2^{b_2}\cdot p_3^{b_3}\cdots p_k^{b_k}$，我们先只考虑 $p_1$（其余可以同理），要使 $\gcd(i, j) = 1$，必然有 $a_1 = 0$ 或 $b_1 = 0$，若 $a_1 = 0$，那么 $b_1$ 有 $x_1 - y_1 + 1$ 种取值；如果 $b_1 = 0$，那么 $a_1$ 有 $y_1 + 1$ 种取值，一共 $x_1 + 1$ 种取值，对其他 $x_i$ 也是如此。
所以满足条件的 $i, j$ 对数为 $\prod(x_i + 1)$，即与约数定理的形式一样。
故原式可以变为
$$\begin{aligned}\sum_{i = 1} ^ n\sum_{j = 1} ^ n\sum_{x | i}\sum_{y | j}[\gcd(x, y) = 1] &= \sum_{i = 1} ^ n\sum_{j = 1} ^ n\sum_{x | i}\sum_{y | j}\sum_{d | x, d | y}\mu(d) \\
&= \sum_{d=1}^n\mu(d)\sum_{x=1}^n[d|x]\sum_{i=1}^n[x|i]\sum\limits_{y=1}^n[d|y]\sum_{j=1}^n[y|j] \\
&= \sum_{d=1}^n\mu(d)(\sum_{x=1}^{\lfloor \frac n d\rfloor}\sum_{i=1}^{\lfloor\frac n d\rfloor}[x|i])^2 \\
&= \sum_{d = 1} ^ n\mu(d)(\sum_{x = 1} ^ {\lfloor \frac n d \rfloor}\lfloor \frac {\lfloor \frac n d \rfloor} {x} \rfloor) ^ 2
\end{aligned}$$
然后 $\mu(d)$ 可以用杜教筛计算，后面的可以分块套分块计算，时间复杂度为 $O(n ^ {\frac 2 3} + n ^ {\frac 3 4})$
代码
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
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
|
#include <bits/stdc++.h>
namespace Task {
#define long long long
const int MAX_SIEVE_N = 1000000;
const int MOD = 1000000007;
int mu[MAX_SIEVE_N + 1], prime[MAX_SIEVE_N], cnt;
bool vis[MAX_SIEVE_N + 1];
inline void fastLinearSieveMethod(const int n) {
mu[1] = 1;
for (register int i = 2; i <= n; i++) {
if (!vis[i]) prime[cnt++] = i, mu[i] = -1;
for (register int j = 0, t; j < n && (t = i * prime[j]) <= n; j++) {
vis[t] = true;
if (i % prime[j] == 0) {
mu[t] = 0;
break;
} else {
mu[t] = -mu[i];
}
}
}
for (register int i = 2; i <= n; i++) mu[i] += mu[i - 1];
}
const int BLOCK_SIZE = 31623, BOUND = 1e9;
int buc1[BLOCK_SIZE + 1], buc2[BLOCK_SIZE + 1];
inline int &get(int x) { return x < BLOCK_SIZE ? buc1[x] : buc2[BOUND / x]; }
inline int sieve(int x) {
if (x <= MAX_SIEVE_N) return mu[x];
register int &cur = get(x);
if (cur) return cur;
register int ret = 1;
for (register int i = 2, pos; i <= x; i = pos + 1) {
pos = x / (x / i);
ret = (ret - (pos - i + 1) * (long)sieve(x / i) % MOD + MOD) % MOD;
}
return cur = ret;
}
inline int getF(int x) {
register int ret = 0;
for (register int i = 1, pos; i <= x; i = pos + 1) {
pos = x / (x / i);
ret = (ret - (pos - i + 1ll) * (x / i) % MOD + MOD) % MOD;
}
return (long)ret * ret % MOD;
}
inline void solve() {
std::ios::sync_with_stdio(false), std::cin.tie(NULL), std::cout.tie(NULL);
register int n, ans = 0;
std::cin >> n;
fastLinearSieveMethod(MAX_SIEVE_N);
for (register int i = 1, pos; i <= n; i = pos + 1) {
pos = n / (n / i),
ans = (ans +
(sieve(pos) - sieve(i - 1) + (long)MOD) % MOD * getF(n / i) %
MOD) %
MOD;
}
std::cout << ans;
}
#undef long
}
int main() {
Task::solve();
return 0;
}
|
