莫反---题解【P1829 [国家集训队]Crash的数字表格 / JZPTAB】
Vocalise
2020-10-21 22:52:59
题解已关。完全自己思路做出来的题目呢...
但可能比较套路（？）
---
求
$$ \sum_{i=1}^n\sum_{j=1}^m\operatorname{lcm}(i,j) $$
其中 $n,m\leq10^7$。
---
直接常规分析：
$$ \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_{d=1}^{\min(n,m)}\sum_{i=1}^{n}\sum_{j=1}^m\frac{i\times j}d[\gcd(i,j)=d] \\
&= \sum_{d=1}^{\min(n,m)}\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}ijd[\gcd(i,j)=1]
\end{aligned} $$
然后没有办法了，反演。
$$ \begin{aligned}
&= \sum_{d=1}^{\min(n,m)}d\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}ij\sum_{t|i,t|j}\mu(t) \\
&= \sum_{d=1}^{\min(n,m)}d\sum_{t=1}^{\min(n,m)/d}\mu(t)\sum_{t|i}^{n/d}i\sum_{t|j}^{m/d}j \\
&= \sum_{d=1}^{\min(n,m)}d\sum_{t=1}^{\min(n,m)/d}\mu(t)t^2\sum_{i=1}^{n/dt}i\sum_{i=1}^{m/dt}j
\end{aligned} $$
记：
$$ f(n) = \sum_{i=1}^ni = \frac {n(n + 1)}2 $$
则
$$ g(n,m) = \sum_{t=1}^{\min(n,m)}\mu(t)t^2f(\lfloor\frac n{t} \rfloor)f(\lfloor\frac mt \rfloor) $$
原式为
$$ \sum_{d=1}^{\min(n,m)}dg(\lfloor\frac nd \rfloor,\lfloor\frac md \rfloor) $$
$f(n)$ 可以 $\Theta(1)$ 求。
而发现 $g(n,m)$ 和答案式都可以数论分块 $\mathcal O(\sqrt n)$ 求，于是...
就可以数论分块中套数论分块了。
考虑计算中的复杂度上界为两次 $\mathcal O(\sqrt n)$ 分块，即 $\mathcal O(n)$。
这部分上界应较松。
不过发现筛 $\mu$ 就是 $\mathcal O(n)$ ，所以总复杂度是 $\mathcal O(n)$ 的没问题。
```cpp
#include <cstdio>
#include <iostream>
#include <cstring>
#include <cmath>
#include <algorithm>
const int N = 10000001;
const int p = 20101009;
inline int read() {
int x = 0, f = 1; char ch = getchar();
while(ch > '9' || ch < '0') { if(ch == '-') f = -1; ch = getchar(); }
do x = x * 10 + ch - 48, ch = getchar(); while(ch >= '0' && ch <= '9');
return x * f;
}
int n,m,mk[N],P[N],t,mu[N],smu[N];
void Seive(int n){
mu[1] = 1;
for(int i = 2;i <= n;i++) {
if(!mk[i]) mu[i] = p - 1, P[++t] = i;
for(int j = 1;j <= t && i * P[j] <= n;j++) {
mk[i * P[j]] = true;
if(i % P[j]) mu[i * P[j]] = p - mu[i];
else mu[i * P[j]] = 0;
}
}
for(int i = 1;i <= n;i++) smu[i] = (smu[i - 1] + 1ll * i * i % p * mu[i] % p) % p;
return;
}
int f(int n) {
return 1ll * n * (n + 1) / 2 % p;
}
int g(int n,int m) {
int ans = 0;
for(int i = 1;i <= n && i <= m;i++) {
int j = std::min(n / (n / i),m / (m / i));
ans = (ans + 1ll * (smu[j] - smu[i - 1] + p) % p * f(n / i) % p * f(m / i) % p) % p;
i = j;
}
return ans;
}
int main() {
n = read(), m = read();
Seive(std::min(n,m));
int ans = 0;
for(int i = 1;i <= n && i <= m;i++) {
int j = std::min(n / (n / i),m / (m / i));
ans = (ans + 1ll * (f(j) - f(i - 1) + p) % p * g(n / i,m / i) % p) % p;
i = j;
}
std::printf("%d\n",ans);
return 0;
}
```