题解:
必要结论:
σ(i∗j)=∑p|i∑q|jp∗jq(gcd(p,q)=1)
证明:
∑p|i∑q|jp∗jq(gcd(p,q)=1)
=∑p|i∑q|jp∗q(gcd(p,jq)=1)
对每一个质因子分类讨论:
设
i=∑puii
,
j=pvii
p=∑plii
,
q=∑prii
1.若
ui>=li>0
,则
ri=vi
,此时可以表示出
p(vi+1)−>(ui+vi)i
。
2.
li=0
,则
vi>=ri>=0
,此时可以表示出
p0−>vii
.
于是可以得出 p∗q 一定可以表示出每个质因子 p0−>(ui+vi)i ,合起来就包括了所有i*j的所有约数。
∴ σ(i∗j)=∑p|i∑q|jp∗jq(gcd(p,q)=1)
有了gcd就可以反演了。
Ans=∑ni=1∑nj=1σ(i∗j)
=∑ni=1∑nj=1∑p|i∑q|jp∗jq(gcd(p,q)=1)
=∑np=1∑nq=1p∗⌊ni⌋∗⌊nj⌋∗(⌊nj⌋+1)/2(gcd(p,q)=1)
=∑np=1∑nq=1q∗⌊ni⌋∗⌊nj⌋∗(⌊nj⌋+1)/2∑d|gcd(p,q)μ(d)
=∑nd=1μ(d)∗d∗∑⌊nd⌋i=1∑⌊nd⌋j=1i∗⌊nd∗i⌋∗⌊nd∗j⌋∗(⌊nd∗j⌋+1)/2
=∑nd=1μ(d)∗d∗(∑⌊nd⌋i=1i∗⌊nd∗i⌋)∗(∑⌊nd⌋i=1⌊nd∗i⌋∗(⌊nd∗i⌋+1)/2)
μ(d)∗d 杜教筛,然后对d分块, ∑⌊nd⌋i=1i∗⌊nd∗i⌋) 和 (∑⌊nd⌋i=1⌊nd∗i⌋∗(⌊nd∗i⌋+1)/2) 继续分块求。
时间复杂度: O(n2/3)
Code:
#include<cmath>
#include<cstdio>
#define ll long long
#define fo(i, x, y) for(ll i = x; i <= y; i ++)
#define min(a, b) ((a) < (b) ? (a) : (b))
using namespace std;
const ll N = 2000000, M = 6516531;
const ll mo = 1e9 + 7, ni_2 = 5e8 + 4;
bool bz[N + 5];
int p[N]; ll mu[N + 5];
ll n, ans;
ll h[M][2];
ll hash(ll x) {
ll y = x % M;
while(h[y][0] != 0 && h[y][0] != x) y = (y + 1) % M;
return y;
}
ll dg(ll n) {
if(n <= N) return mu[n];
ll y = hash(n);
if(h[y][0] == n) return h[y][1];
h[y][0] = n;
ll s = 1;
fo(i, 1, n) {
ll j = n / (n / i);
s -= dg(n / i) * ((j - i + 1) * (i + j) / 2 % mo) % mo;
i = j;
}
s = (s % mo + mo) % mo;
h[y][1] = s;
return s;
}
ll Solve(ll n) {
ll p = 0, q = 0;
fo(i, 1, n) {
ll j = n / (n / i);
p += (j - i + 1) * (i + j) / 2 % mo * (n / i);
q += (n / i) * (n / i + 1) / 2 % mo * (j - i + 1) % mo;
i = j;
}
return p % mo * (q % mo) % mo;
}
int main() {
mu[1] = 1;
fo(i, 2, N) {
if(!bz[i]) p[++ p[0]] = i, mu[i] = -1;
fo(j, 1, p[0]) {
ll k = i * p[j];
if(k > N) break;
bz[k] = 1;
if(i % p[j] == 0) {
mu[k] = 0; break;
}
mu[k] = -mu[i];
}
}
fo(i, 1, N) mu[i] = (mu[i] * i + mu[i - 1] + mo) % mo;
scanf("%lld", &n);
fo(d, 1, n) {
ll d_ = n / (n / d);
ans += Solve(n / d) * (dg(d_) - dg(d - 1) + mo) % mo;
d = d_;
}
printf("%lld", ans % mo);
}