题目描述

给定$n$,求$\sum_{i=1}^{n}\sum_{j=1}^{n}\sigma_1(ij)$

其中$2 \le n \le 10^{9}$

题解

首先有

$$
\sigma_1(ij)=\sum_{x \mid i}\sum_{y \mid j}[(x,y)=1]\frac{xj}{y}
$$

那么可以得到

$$
\begin{align}
\sum_{i=1}^{n}\sum_{j=1}^{n}\sigma_1(ij)
&=\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{x \mid i}\sum_{y \mid j}[(x,y)=1]\frac{xj}{y} \\
&=\sum_{d=1}^{n}\mu(d)\sum_{i=1}^{n}\sum_{j=1}^{n}\sum_{x \mid i}\sum_{y \mid j}[d \mid (x,y)]\frac{xj}{y} \\
&=\sum_{d=1}^{n}\mu(d)\sum_{d \mid x}\sum_{d \mid y}\frac{x}{y}\sum_{i=1 \wedge x \mid i}^{n} \sum_{j=1 \wedge y \mid j}^{n}j \\
&=\sum_{d=1}^{n}\mu(d)\sum_{d \mid x}\sum_{d \mid y}\frac{x}{y}\lfloor \frac{n}{x} \rfloor y \frac{(1 + \lfloor \frac{n}{y} \rfloor)\lfloor \frac{n}{y} \rfloor}{2} \\
&=\sum_{d=1}^{n}\mu(d) (\sum_{d \mid x}x\lfloor \frac{n}{x} \rfloor) (\sum_{d \mid y}\frac{(1 + \lfloor \frac{n}{y} \rfloor)\lfloor \frac{n}{y} \rfloor}{2}) \\
&=\sum_{d=1}^{n}\mu(d) (\sum_{x=1}^{\lfloor \frac{n}{d} \rfloor}xd \lfloor \frac{n}{dx} \rfloor) (\sum_{y=1}^{\lfloor \frac{n}{d} \rfloor}\frac{(1+\lfloor \frac{n}{dy} \rfloor)\lfloor \frac{n}{dy} \rfloor}{2}) \\
&=\sum_{d=1}^{n}\mu(d)d (\sum_{x=1}^{\lfloor \frac{n}{d} \rfloor}x \lfloor \frac{n}{dx} \rfloor) (\sum_{y=1}^{\lfloor \frac{n}{d} \rfloor}\frac{(1+\lfloor \frac{n}{dy} \rfloor)\lfloor \frac{n}{dy} \rfloor}{2}) \\
\end{align}
$$

设$p(n)=\sum_{x=1}^{n} x\lfloor \frac{n}{x} \rfloor$,相当于每个$x$被其所有$[1,n]$之间的倍数各做了一次贡献,即$p(n)=\sum_{x=1}^{n}\sigma_1(x)$

设$q(n)=\sum_{x=1}^{n}\frac{(1+\lfloor \frac{n}{x} \rfloor)\lfloor \frac{n}{x} \rfloor}{2}$

本质上是

$$
\begin{align}
q(n)&=\sum_{x=1}^{n}\frac{(1+\lfloor \frac{n}{x} \rfloor)\lfloor \frac{n}{x} \rfloor}{2} \\
&=\sum_{x=1}^{n}\sum_{j=1}^{\lfloor \frac{n}{x} \rfloor}j \\
&=\sum_{j=1}^{n}j\sum_{x=1}^{\lfloor \frac{n}{j} \rfloor} 1 \\
&=\sum_{x=1}^{n}x\lfloor \frac{n}{x} \rfloor \\
&=\sum_{x=1}^{n}\sigma_1(x)
\end{align}
$$

所以

$$
\begin{align}
\sum_{i=1}^{n}\sum_{j=1}^{n} \sigma_1(ij)
&= \sum_{d=1}^{n}\mu(d)d (\sum_{x=1}^{\lfloor \frac{n}{d} \rfloor}\sigma_1(x))^2 \\
&= \sum_{d=1}^{n}\mu(d)d \times h(\lfloor \frac{n}{d} \rfloor)^2 \\
\end{align}
$$

不妨设$f(d)=\mu(d)d$,则有$(f \times id)(n) = \sum_{d \mid n} \mu(d)d \frac{n}{d}=n \epsilon(n)=\epsilon(n)$

即$F(n)=1-\sum_{d=2}^{n}dF(\lfloor \frac{n}{d} \rfloor)$

对于$f(d)$还需要求出前$n^{\frac{2}{3}}$项

首先有$f(1)=\mu(1)1=1$

设当前枚举到了$i$和素数$p_j$

如果$i$是素数,则$f(i)=\mu(i)i=-i$

如果$p_j \not \mid i$,则$f(i)=\mu(i \times p_j)i \times p_j=-\mu(i)i \times p_j=-f(i) \times p_j$

如果$p_j \mid i$,则$f(i)=\mu(i \times p_j)i \times p_j=\mu(i’ \times p_j^2)i \times p_j=0$

对于$h(\lfloor \frac{n}{d} \rfloor)$,可以线性筛预处理前$n^{\frac{2}{3}}$项

设$e(n)$表示$n$的最小素因子的幂次和,$d(n)$表示$n$的约数和(即$\sigma_1(n)$)

设$n=\Pi p_i^{k_i}$,则有$d(n)=\Pi_{i} \sum_{j=0}^{k_i}p_i^{j}$ ,$e(n)=\sum_{j=0}^{k_{min}}p_{min}^{j}$

设当前枚举到了$i$和素数$p_j$

若$i$是素数,则有$d(i)=1+i$,$e(i)=1+i$

若$p_j \not \mid i$,此时$i \times p_j$的最小素因子是$p_j$,即$d(i \times p_j)=d(i) \times (1+p_j)$,$e(i \times p_j)=1+p_j$

若$p_j \mid i$,此时$i \times p_j$和$i$的最小素因子都是$p_j$,则有$d(i \times p_j)=\frac{d(i)}{e(i)}(p_j \times e(i)+1)$,$e(i \times p_j)=p_j \times e(i)+1$

对于后$n^{\frac{1}{3}}$项,可以根据$h(n)=\sum_{x=1}^{n}x \lfloor \frac{n}{x} \rfloor$进行数论分块,$O(\sqrt n)$计算

最终计算$h(n)$的时间复杂度为$O(n^{\frac{2}{3}})$

由于$\lfloor \frac{n}{d} \rfloor$的取值只有$\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
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
85
86
87
88
#include "bits/stdc++.h"
using namespace std;
typedef long long ll;

const int N = 1e6, mod = 1e9 + 7;

ll e[N + 10], d[N + 10], muid[N + 10], mu[N + 10];

int pri[N + 10], tot, vis[N + 10];

ll n, ans;

map<ll, ll> val_sig;

ll sum_sig1(ll n) {
if(n <= N) {
return d[n];
} else {
if(val_sig.find(n) != val_sig.end()) {
return val_sig[n];
} else {
ll res = 0;
for(ll i = 1, j ; i <= n ; i = j + 1) {
j = n / (n / i);
res = (res + (i + j) * (j - i + 1) / 2 * (n / i) % mod) % mod;
}
return val_sig[n] = res;
}
}
}

map<ll, ll> val_muid;

ll sum_muid(ll n) {
if(n <= N) {
return muid[n];
} else {
if(val_muid.find(n) != val_muid.end()) {
return val_muid[n];
} else {
ll res = 1;
for(ll i = 2, j ; i <= n ; i = j + 1) {
j = n / (n / i);
res = (res - (i + j) * (j - i + 1) / 2 % mod * sum_muid(n / i)) % mod;
}
return val_muid[n] = res;
}
}
}

int main() {
scanf("%lld", &n);

d[1] = 1;
mu[1] = 1;
muid[1] = 1;
for(int i = 2 ; i <= N ; ++ i) {
if(!vis[i]) {
d[i] = 1 + i;
e[i] = 1 + i;
mu[i] = -1;
pri[++ tot] = i;
}
for(int j = 1 ; j <= tot && (ll) i * pri[j] <= N ; ++ j) {
vis[i * pri[j]] = 1;
if(i % pri[j] == 0) {
d[i * pri[j]] = d[i] / e[i] * (pri[j] * e[i] + 1);
e[i * pri[j]] = pri[j] * e[i] + 1;
break;
} else {
d[i * pri[j]] = d[i] * (1 + pri[j]);
e[i * pri[j]] = 1 + pri[j];
mu[i * pri[j]] = -mu[i];
}
}
d[i] = (d[i] + d[i - 1]) % mod;
muid[i] = (muid[i - 1] + mu[i] * i) % mod;
}
for(ll i = 1, j ; i <= n ; i = j + 1) {
j = n / (n / i);
ll sig = sum_sig1(n / i) % mod;
sig = sig * sig % mod;
ll mui = (sum_muid(j) - sum_muid(i - 1)) % mod;
ans = (ans + mui % mod * sig % mod) % mod;
}

printf("%lld\n", (ans % mod + mod) % mod);
}