「国家集训队」Crash的数字表格
题目描述
求(对 \(20101009\) 取模，\(n,m\le10^7\) )
\[\sum_{i=1}^n\sum_{j=1}^m\operatorname{lcm}(i,j)
\]
大体思路
推式子:
\[\begin{aligned}\sum_{i=1}^n\sum_{j=1}^m\operatorname{lcm}(i,j)
&=\sum_{i=1}^n\sum_{j=1}^m\frac{i\times j}{\gcd(i,j)} \\
&=\sum_{i=1}^n\sum_{j=1}^m\sum_{d|i,d|j,\gcd(i/d,j/d)=1}\frac{i\times j}{d} \\
&=\sum_{d=1}^{\min(n,m)}\times d\times\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}[\gcd(i,j)=1]\times i\times j\end{aligned}\]
把式子后面那一大堆设为 \(sum(n,m)\) ：
\[sum(n,m)=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=1]\times i\times j
\]
考虑化简一下 \(sum\) :
\[\begin{aligned}sum(n,m) &= \sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=1]\times i\times j\\
&=\sum_{i=1}^n\sum_{j=1}^m\sum_{d|\gcd(i,j)}\mu(d)\times i\times j \\
&=\sum_{d=1}^{\min(n,m)}\mu(d)\times d^2\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}i\times j\end{aligned}\]
可以发现 \(sum\) 后面那一大堆(设为 \(g(n,m)\) )可以 \(O(1)\) 求：
\[\begin{aligned}g(n,m)&=\sum_{i=1}^n\sum_{j=1}^m i\times j \\
&=\frac{n\times(n+1)}{2}\times \frac{m\times(m+1)}{2}\end{aligned}\]
那么 \(sum(n,m)\) 可以化为：
\[sum(n,m)=\sum_{d=1}^{\min(n,m)}\mu(d)\times d^2\times g(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)
\]
这个可以数论分块 \(\lfloor\frac{n}{\lfloor\frac{n}{d}\rfloor}\rfloor\) 求。
再回到定义 \(sum\) 的地方，那么：
\[Ans=\sum_{d=1}^{\min(n,m)}\times d\times sum(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor)
\]
好像这个还是可以数论分块 \(QwQ\)
至此这道题就解决了。
细节注意事项
- \(long\ long\)一定要开呀。
- 不要写挂呀！！！
参考代码
/*--------------------------------
Code name: crash.cpp
Author: The Ace Bee
This code is made by The Ace Bee
--------------------------------*/
#include <cstdio>
#define rg register
#define int long long
#define fileopen(x) \
freopen(x".in", "r", stdin); \
freopen(x".out", "w", stdout);
#define fileclose \
fclose(stdin); \
fclose(stdout);
const int mod = 20101009;
const int MAXN = 10000010;
inline int min(int a, int b) { return a < b ? a : b; }
inline int read() {
int s = 0; bool f = false; char c = getchar();
while (c < '0' || c > '9') f |= (c == '-'), c = getchar();
while (c >= '0' && c <= '9') s = (s << 3) + (s << 1) + (c ^ 48), c = getchar();
return f ? -s : s;
}
int vis[MAXN], mu[MAXN];
int num, pri[MAXN], sum[MAXN];
inline void seive() {
mu[1] = 1;
for (rg int i = 2; i < MAXN; ++i) {
if (!vis[i]) mu[i] = -1, pri[++num] = i;
for (rg int j = 1; j <= num && i * pri[j] < MAXN; ++j) {
vis[i * pri[j]] = 1;
if (i % pri[j]) mu[i * pri[j]] = - mu[i];
else { mu[i * pri[j]] = 0; break; }
}
}
for (rg int i = 1; i < MAXN; ++i)
sum[i] = (sum[i - 1] + 1ll * i * i % mod * (mu[i] + mod) % mod) % mod;
}
inline int g(int n, int m)
{ return 1ll * n * (n + 1) / 2 % mod * (m * (m + 1) / 2 % mod) % mod; }
inline int f(int n, int m) {
int res = 0;
for (rg int i = 1, j; i <= min(n, m); i = j + 1) {
j = min(n / (n / i), m / (m / i));
res = (res + 1ll * (sum[j] - sum[i - 1] + mod) * g(n / i, m / i) % mod) % mod;
}
return res;
}
inline int solve(int n, int m) {
int res = 0;
for (rg int i = 1, j; i <= min(n, m); i = j + 1) {
j = min(n / (n / i), m / (m / i));
res = (res + 1ll * (j - i + 1) * (i + j) / 2 % mod * f(n / i, m / i) % mod) % mod;
}
return res;
}
signed main() {
// fileopen("crash");
seive();
int n = read(), m = read();
printf("%lld\n", solve(n, m));
// fileclose;
return 0;
}
完结撒花\(qwq\)
