∑
i
=
1
n
∑
j
=
1
n
i
∗
j
∗
g
c
d
(
i
,
j
)
\sum_{i = 1}^n\sum_{j = 1}^ni * j * gcd(i,j)
i=1∑nj=1∑ni∗j∗gcd(i,j)
=
∑
d
=
1
n
∑
i
=
1
⌊
n
d
⌋
∑
j
=
1
⌊
n
d
⌋
i
∗
j
∗
d
3
∗
[
g
c
d
(
i
,
j
)
=
1
]
=\sum_{d = 1}^{n}\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j = 1}^{\lfloor\frac{n}{d}\rfloor} i * j * d^3 * [gcd(i,j) = 1]
=d=1∑ni=1∑⌊dn⌋j=1∑⌊dn⌋i∗j∗d3∗[gcd(i,j)=1]
=
∑
d
=
1
n
∑
i
=
1
⌊
n
d
⌋
∑
j
=
1
⌊
n
d
⌋
i
∗
j
∗
d
3
∗
∑
p
∣
g
c
d
(
i
,
j
)
μ
(
p
)
=\sum_{d = 1}^{n}\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j = 1}^{\lfloor\frac{n}{d}\rfloor} i * j * d^3 * \sum_{p | gcd(i,j)}\mu(p)
=d=1∑ni=1∑⌊dn⌋j=1∑⌊dn⌋i∗j∗d3∗p∣gcd(i,j)∑μ(p)
=
∑
d
=
1
n
∑
p
=
1
⌊
n
d
⌋
μ
(
p
)
∑
i
=
1
⌊
n
d
∗
p
⌋
∑
j
=
1
⌊
n
d
∗
p
⌋
i
∗
j
∗
d
3
∗
p
2
=\sum_{d = 1}^{n}\sum_{p =1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\sum_{i = 1}^{\lfloor\frac{n}{d * p}\rfloor}\sum_{j = 1}^{\lfloor\frac{n}{d * p}\rfloor} i * j * d^3*p^2
=d=1∑np=1∑⌊dn⌋μ(p)i=1∑⌊d∗pn⌋j=1∑⌊d∗pn⌋i∗j∗d3∗p2
=
∑
T
=
1
n
∑
d
∣
T
μ
(
T
d
)
∗
d
3
∗
(
T
d
)
2
∗
f
(
⌊
n
T
⌋
)
2
,
f
(
n
)
=
∑
i
=
1
n
i
=\sum_{T = 1}^{n}\sum_{d | T}\mu(\frac{T}{d})*d^3*(\frac{T}{d})^2*f({\lfloor\frac{n}{T}\rfloor})^2,f(n) = \sum_{i = 1}^ni
=T=1∑nd∣T∑μ(dT)∗d3∗(dT)2∗f(⌊Tn⌋)2,f(n)=i=1∑ni
=
∑
T
=
1
n
T
2
∑
d
∣
T
μ
(
T
d
)
∗
d
∗
f
(
⌊
n
T
⌋
)
2
=\sum_{T = 1}^{n}T^2\sum_{d | T}\mu(\frac{T}{d})*d*f({\lfloor\frac{n}{T}\rfloor})^2
=T=1∑nT2d∣T∑μ(dT)∗d∗f(⌊Tn⌋)2
=
∑
T
=
1
n
T
2
ϕ
(
T
)
∗
f
(
⌊
n
T
⌋
)
2
=\sum_{T = 1}^{n}T^2\phi(T)*f({\lfloor\frac{n}{T}\rfloor})^2
=T=1∑nT2ϕ(T)∗f(⌊Tn⌋)2
令
f
(
T
)
=
T
2
∗
ϕ
(
T
)
f(T) = T ^2 * \phi(T)
f(T)=T2∗ϕ(T)
n
<
=
1
0
7
n <= 10^7
n<=107的话可以用线性筛预处理
f
f
f函数前缀和,然后分块。
n 出到了
1
0
11
10^{11}
1011可以用杜教筛筛出
f
f
f函数的前缀和,然后分块
杜教筛的大致过程:选一个 积性函数
g
g
g,构造
h
=
g
∗
f
h = g * f
h=g∗f(
∗
*
∗代表狄利克雷卷积),可知 h 也是一个积性函数。
h
(
x
)
=
∑
d
∣
x
g
(
d
)
∗
f
(
x
d
)
h(x) = \sum_{d | x} g(d) * f(\frac{x}{d})
h(x)=d∣x∑g(d)∗f(dx)
∑
i
=
1
n
h
(
i
)
=
∑
i
=
1
n
∑
d
∣
i
g
(
d
)
∗
f
(
i
d
)
\sum_{i = 1}^nh(i) = \sum_{i = 1}^n\sum_{d | i} g(d) * f(\frac{i}{d})
i=1∑nh(i)=i=1∑nd∣i∑g(d)∗f(di)
=
∑
d
=
1
n
∑
i
=
1
⌊
n
d
⌋
g
(
d
)
∗
f
(
i
)
=\sum_{d = 1}^n\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor} g(d) * f(i)
=d=1∑ni=1∑⌊dn⌋g(d)∗f(i)
令
S
(
n
)
=
∑
i
=
1
n
f
(
i
)
:
∑
i
=
1
n
h
(
i
)
=
∑
d
=
1
n
∑
i
=
1
⌊
n
d
⌋
g
(
d
)
∗
f
(
i
)
=
∑
d
=
1
n
g
(
d
)
∗
S
(
⌊
n
d
⌋
)
=
g
(
1
)
∗
S
(
n
)
+
∑
d
=
2
n
g
(
d
)
∗
S
(
⌊
n
d
⌋
)
令S(n) = \sum_{i = 1}^nf(i):\sum_{i = 1}^nh(i) =\sum_{d = 1}^n\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor} g(d) * f(i) =\sum_{d = 1}^ng(d)*S({\lfloor\frac{n}{d}\rfloor}) = g(1) * S(n) + \sum_{d = 2}^ng(d)*S({\lfloor\frac{n}{d}\rfloor})
令S(n)=i=1∑nf(i):i=1∑nh(i)=d=1∑ni=1∑⌊dn⌋g(d)∗f(i)=d=1∑ng(d)∗S(⌊dn⌋)=g(1)∗S(n)+d=2∑ng(d)∗S(⌊dn⌋)
移
项
:
g
(
1
)
∗
S
(
n
)
=
∑
i
=
1
n
h
(
i
)
−
∑
d
=
2
n
g
(
d
)
∗
S
(
⌊
n
d
⌋
)
移项:g(1) * S(n) = \sum_{i = 1}^nh(i) - \sum_{d = 2}^ng(d)*S({\lfloor\frac{n}{d}\rfloor})
移项:g(1)∗S(n)=i=1∑nh(i)−d=2∑ng(d)∗S(⌊dn⌋)若
h
(
i
)
h(i)
h(i)的前缀和可以在O(1)时间求得,则
S
(
n
)
S(n)
S(n)可以通过分块加记忆化搜索来求,因此要选一个积性函数
f
f
f使得
h
=
g
∗
f
h = g * f
h=g∗f且
h
h
h 的前缀和容易求。
观
察
g
∗
f
:
∑
d
∣
x
g
(
d
)
∗
f
(
x
d
)
=
∑
d
∣
x
g
(
d
)
∗
(
x
d
)
2
∗
ϕ
(
x
d
)
观察g * f:\sum_{d | x}g(d) * f(\frac{x}{d}) = \sum_{d | x}g(d) * ({\frac{x}{d}})^2*\phi(\frac{x}{d})
观察g∗f:d∣x∑g(d)∗f(dx)=d∣x∑g(d)∗(dx)2∗ϕ(dx)
当
g
(
i
)
g(i)
g(i) 取
i
2
i^2
i2时(方便去掉分母里的
d
2
d^2
d2)
g
∗
f
=
∑
d
∣
x
d
2
∗
ϕ
(
x
d
)
=
x
3
,
前
缀
和
公
式
为
:
x
2
∗
(
x
+
1
)
2
4
g * f = \sum_{d | x}d ^ 2 * \phi(\frac{x}{d}) =x ^ 3,前缀和公式为:\frac{x^2*(x + 1)^2}{4}
g∗f=d∣x∑d2∗ϕ(dx)=x3,前缀和公式为:4x2∗(x+1)2
∑
i
=
1
n
i
2
=
i
∗
(
i
+
1
)
∗
(
2
∗
i
+
1
)
6
\sum_{i = 1}^n i ^2 = \frac{i * (i + 1) * (2*i + 1)}{6}
i=1∑ni2=6i∗(i+1)∗(2∗i+1)
于是就可以杜教筛了,注意记忆化,先用线性筛预处理
n
\sqrt n
n的前缀和,记忆化可以用map 也可以 unordered_map 或自己的hash,由于洛谷unordered_map编译不通过,这里还是用map
代码:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 6e6;
map<ll,ll> mp;
bool ispri[maxn + 100];
int pri[maxn + 100],phi[maxn + 100],mu[maxn + 100];
ll sum[maxn + 100];
ll n,mod;
ll inv2,inv4,inv6;
ll fpow(ll a,ll b) {
ll r = 1;
while(b) {
if(b & 1) r = r * a % mod;
a = a * a % mod;
b >>= 1;
}
return r;
}
void sieve(int n) {
ispri[0] = ispri[1] = false;
pri[0] = 0;mu[1] = 1;phi[1] = 1;
for(int i = 2; i <= n; i++) {
if(!ispri[i]) pri[++pri[0]] = i,mu[i] = -1,phi[i] = i - 1;
for(int j = 1; j <= pri[0] && i * pri[j] <= n; j++) {
ispri[i * pri[j]] = true;
if(i % pri[j] == 0) {
phi[i * pri[j]] = phi[i] * pri[j];
break;
}
phi[i * pri[j]] = phi[i] * (pri[j] - 1);
mu[i * pri[j]] = -mu[i];
}
}
sum[0] = 0;
for(int i = 1; i <= n; i++)
sum[i] = (sum[i - 1] + phi[i] * (1ll * i * i % mod) % mod) % mod;
}
ll cal0(ll x) {
x %= mod;
return x * (x + 1) % mod * inv2 % mod;
}
ll cal1(ll x) {
x %= mod;
return x * x % mod * (x + 1) % mod * (x + 1) % mod * inv4 % mod;
}
ll cal2(ll x) {
x %= mod;
return x * (x + 1) % mod * (2 * x + 1) % mod * inv6 % mod;
}
ll getsum(ll p) {
if(p <= maxn) return sum[p];
if(mp[p]) return mp[p];
ll l,r;
ll res = cal1(p);
for(l = 2; l <= p; l = r + 1) {
r = (p / (p / l));
ll t = getsum(p / l);
res = (res + mod - t * ((cal2(r) - cal2(l - 1) + mod) % mod) % mod) % mod;
}
return mp[p] = (res + mod) % mod;
}
int main() {
scanf("%lld%lld",&mod,&n);
sieve(6000000);
inv2 = fpow(2,mod - 2);
inv4 = fpow(4,mod - 2);
inv6 = fpow(6,mod - 2);
ll l,r;
ll ans = 0;
for(l = 1; l <= n; l = r + 1) {
r = (n / (n / l));
ll t = (getsum(r) - getsum(l - 1) + mod) % mod;
ll x = cal1(n / l);
ans = (ans + x * t % mod) % mod;
}
printf("%lld\n",(ans + mod) % mod);
return 0;
}