洛谷 P1829 [国家集训队]Crash的数字表格 / JZPTAB
题意
求
∑
i
=
1
n
∑
j
=
1
m
lcm
(
i
,
j
)
(
m
o
d
20101009
)
\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\text{lcm}(i,j)(\bmod 20101009)
i=1∑nj=1∑mlcm(i,j)(mod20101009)
思路
容易想到原式等价于
∑
i
=
1
n
∑
j
=
1
m
i
∗
j
gcd
(
i
,
j
)
\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{i* j}{\gcd(i,j)}
i=1∑nj=1∑mgcd(i,j)i∗j
枚举
i
,
j
i,j
i,j的最大公约数
d
d
d,显然
gcd
(
i
d
,
j
d
)
=
1
\gcd(\frac id,\frac jd)=1
gcd(di,dj)=1,即
i
d
\frac id
di和
j
d
\frac jd
dj互质
∑
i
=
1
n
∑
j
=
1
m
∑
d
∣
i
,
d
∣
j
,
gcd
(
i
d
,
j
d
)
=
1
i
∗
j
d
\sum\limits_{i=1}^{n}\sum\limits_{j=1}^m\sum\limits_{d|i,d|j,\gcd(\frac id,\frac jd)=1}\frac{i*j}d
i=1∑nj=1∑md∣i,d∣j,gcd(di,dj)=1∑di∗j
变换求和顺序
∑
d
=
1
n
d
∑
i
=
1
⌊
n
d
⌋
∑
j
=
1
⌊
m
d
⌋
[
gcd
(
i
,
j
)
=
1
]
i
∗
j
\sum\limits_{d=1}^{n}d\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}[\gcd(i,j)=1]i*j
d=1∑ndi=1∑⌊dn⌋j=1∑⌊dm⌋[gcd(i,j)=1]i∗j
记
s
u
m
(
n
,
m
)
=
∑
i
=
1
n
∑
j
=
1
m
[
gcd
(
i
,
j
)
=
1
]
i
∗
j
sum(n,m)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)=1]i*j
sum(n,m)=i=1∑nj=1∑m[gcd(i,j)=1]i∗j
对其进行化简,用
ε
(
gcd
(
i
,
j
)
)
\varepsilon(\gcd(i,j))
ε(gcd(i,j))替换
[
gcd
(
i
,
j
)
=
1
]
[\gcd(i,j)=1]
[gcd(i,j)=1]
∑
i
=
1
n
∑
j
=
1
m
∑
d
∣
gcd
(
i
,
j
)
μ
(
d
)
∗
i
∗
j
\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{d|\gcd(i,j)}\mu(d)*i*j
i=1∑nj=1∑md∣gcd(i,j)∑μ(d)∗i∗j
转化为首先枚举约数
∑
d
=
1
min
(
n
,
m
)
∑
d
∣
i
n
∑
d
∣
j
m
μ
(
d
)
∗
i
∗
j
\sum\limits_{d=1}^{\min(n,m)}\sum\limits_{d|i}^{n}\sum\limits_{d|j}^{m}\mu(d)*i*j
d=1∑min(n,m)d∣i∑nd∣j∑mμ(d)∗i∗j
设
i
=
i
′
∗
d
,
j
=
j
′
∗
d
i=i'*d,j=j'*d
i=i′∗d,j=j′∗d,则可以进一步转化
∑
d
=
1
min
(
n
,
m
)
μ
(
d
)
∗
d
2
∗
∑
i
=
1
⌊
n
d
⌋
∑
j
=
1
⌊
m
d
⌋
i
∗
j
\sum\limits_{d=1}^{\min(n,m)}\mu(d)*d^2*\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i*j
d=1∑min(n,m)μ(d)∗d2∗i=1∑⌊dn⌋j=1∑⌊dm⌋i∗j
前半段可以处理前缀和,后半段可以
O
(
1
)
O(1)
O(1)求,设
Q
(
n
,
m
)
=
∑
i
=
1
n
∑
j
=
1
m
i
∗
j
=
n
∗
(
n
+
1
)
2
∗
m
∗
(
m
+
1
)
2
Q(n,m)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}i*j=\frac{n*(n+1)}{2}*\frac{m*(m+1)}{2}
Q(n,m)=i=1∑nj=1∑mi∗j=2n∗(n+1)∗2m∗(m+1)
显然可以
O
(
1
)
O(1)
O(1)求解
到现在
s
u
m
(
n
,
m
)
=
∑
d
=
1
min
(
n
,
m
)
μ
(
d
)
∗
d
2
∗
Q
(
⌊
n
d
⌋
,
⌊
m
d
⌋
)
sum(n,m)=\sum\limits_{d=1}^{\min(n,m)}\mu(d)*d^2*Q(\lfloor\frac nd \rfloor,\lfloor\frac md\rfloor)
sum(n,m)=d=1∑min(n,m)μ(d)∗d2∗Q(⌊dn⌋,⌊dm⌋)
可以用数论分块求解
回带到原式中
∑
d
=
1
min
(
n
,
m
)
d
∗
s
u
m
(
⌊
n
d
⌋
,
⌊
m
d
⌋
)
\sum\limits_{d=1}^{\min(n, m)}d*sum(\lfloor\frac nd \rfloor,\lfloor\frac md\rfloor)
d=1∑min(n,m)d∗sum(⌊dn⌋,⌊dm⌋)
又可以数论分块求解了
然后就做完啦
代码
/*
Author:loceaner
*/
#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define int long long
using namespace std;
const int A = 1e7 + 11;
const int B = 1e6 + 11;
const int mod = 20101009;
const int inf = 0x3f3f3f3f;
inline int read() {
char c = getchar(); int x = 0, f = 1;
for( ; !isdigit(c); c = getchar()) if(c == '-') f = -1;
for( ; isdigit(c); c = getchar()) x = x * 10 + (c ^ 48);
return x * f;
}
bool vis[A];
int n, m, mu[A], p[B], sum[A], cnt;
void getmu() {
mu[1] = 1;
int k = min(n, m);
for (int i = 2; i <= k; i++) {
if (!vis[i]) p[++cnt] = i, mu[i] = -1;
for (int j = 1; j <= cnt && i * p[j] <= k; ++j) {
vis[i * p[j]] = 1;
if (i % p[j] == 0) break;
mu[i * p[j]] = -mu[i];
}
}
for (int i = 1; i <= k; i++) sum[i] = (sum[i - 1] + i * i % mod * mu[i]) % mod;
}
int Sum(int x, int y) {
return (x * (x + 1) / 2 % mod) * (y * (y + 1) / 2 % mod) % mod;
}
int solve2(int x, int y) {
int res = 0;
for (int i = 1, j; i <= min(x, y); i = j + 1) {
j = min(x / (x / i), y / (y / i));
res = (res + 1LL * (sum[j] - sum[i - 1] + mod) * Sum(x / i, y / i) % mod) % mod;
}
return res;
}
int solve(int x, int y) {
int res = 0;
for (int i = 1, j; i <= min(x, y); i = j + 1) {
j = min(x / (x / i), y / (y / i));
res = (res + 1LL * (j - i + 1) * (i + j) / 2 % mod * solve2(x / i, y / i) % mod) % mod;
}
return res;
}
signed main() {
n = read(), m = read();
getmu();
cout << solve(n, m) << '\n';
}