Description
求 $$ \sum\limits_{i=1}^n\sum\limits_{j=1}^m\operatorname{lcm}(i,j) $$ $n,m\leq10^7$。
Solution
这道题算是一道比较有代表性的题目了，综合了各种 trick。
trick1：lcm 转化为 gcd： $$ \sum\limits_{i=1}^n\sum\limits_{j=1}^m\operatorname{lcm}(i,j)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{i\times j}{\gcd(i,j)} $$ trick2：枚举分母： $$ ans=\sum\limits_{d=1}^{\min(n,m)}\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{i\times j}{d}[\gcd(i,j)=d] $$ trick3：把 $[\gcd(i,j)=d]$ 转化为 $[\gcd(\dfrac{i}{d},\dfrac{j}{d})=1]$，同时改变枚举的 $i,j$（$i$ 变为 $i/d$，$j$ 同理）。 $$ \begin{aligned} ans&=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{i\times j}{d\times d}[\gcd(\frac{i}{d},\frac{j}{d})=1]\\ &=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{i=1}^{\lfloor \frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j[\gcd(i,j)=1] \end{aligned} $$ trick4：$[n=1]$ 用一个结论转化为莫比乌斯函数的形式：$\sum_{d\mid n}\mu(d)=[n=1]$（证明不讲了）。 $$ ans=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\sum\limits_{k\mid\gcd(i,j)}\mu(k) $$ trick5：优先枚举因数，然后把枚举因数转化成枚举倍数放到后面的限制条件里。 $$ \begin{aligned} ans&=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{k\mid\gcd(i,j)}\mu(k)\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\\ &=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{k=1}^{\lfloor\min(n,m)/d\rfloor}\mu(k)\sum\limits_{i=1,k\mid i}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1,k\mid j}^{\lfloor\frac{m}{d}\rfloor}i\times j \end{aligned} $$ trick6：$1\sim n$ 中 $x$ 的倍数有 $\lfloor\dfrac{n}{x}\rfloor$ 个，那么有 $\sum_{i=1,x\mid i}^nt=\sum_{i=1}^{\lfloor n/x\rfloor}tx$（不一定完全一样，但可以用这个方法化简）。 $$ \begin{aligned} ans&=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{k=1}^{\lfloor\min(n,m)/d\rfloor}\mu(k)\times k^2\sum\limits_{i=1,k\mid i}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1,k\mid j}^{\lfloor\frac{m}{d}\rfloor}\frac{i\times j}{k\times k}\\ &=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{k=1}^{\lfloor\min(n,m)/d\rfloor}\mu(k)\times k^2\sum\limits_{i=1}^{\lfloor\frac{n}{dk}\rfloor}i\sum\limits_{j=1}^{\lfloor\frac{m}{dk}\rfloor}j \end{aligned} $$ 到这一步，$k^2$ 后面的两个东西是等差数列求和的形式，可以 $O(1)$ 得出，然后看到除法+下取整，而且 $\mu(k)\times k^2$ 的前缀和可以预处理得出，那么这玩意就是个数论分块：
设 $$ sum(n)=\sum\limits_{i=1}^ni $$ $$ f(x,y)=\sum\limits_{k=1}^{\lfloor\min(n,m)\rfloor}\mu(k)\times k^2\times sum(\lfloor\frac{n}{k}\rfloor)\times sum(\lfloor\frac{m}{k}\rfloor) $$ 那么 $$ ans=\sum\limits_{d=1}^{\min(n,m)}d\times f(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor) $$ 诶，这又是个数论分块的形式？？？
那这就做完了。
时间复杂度似乎应该是 $O(n^{0.75}+n)$（后面是预处理的复杂度）。
不过我们的式子是可以接着往下推的，还有 trick 可以讲。
可以发现我们上面的式子中有 $k,d,dk$ 三个未知数，我们可以枚举其中一个来消去一个未知数。
trick7：枚举一个未知数来达到减小未知数的目的（这个好像比较灵活，我讲不太明白）。
设 $s=dk$。 $$ ans=\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{s=1,d\mid s}^{\min(n,m)}\mu(\frac{s}{d})\times(\frac{s}{d})^2\times sum(\frac{n}{s})\times sum(\frac{m}{s}) $$ trick8：这个算是很通用的 trick 了，换枚举顺序，同时换枚举的限制，让式子更舒服（？）。 $$ \begin{aligned} ans&=\sum\limits_{s=1}^{\min(n,m)}sum(\frac{n}{s})\times sum(\frac{m}{s})\sum\limits_{d=1,d\mid s}d\times\mu(\frac{s}{d})\times(\frac{s}{d})^2\\ &=\sum\limits_{s=1}^{\min(n,m)}s\times sum(\frac{n}{s})\times sum(\frac{m}{s})\sum\limits_{d=1,d\mid s}\mu(\frac{s}{d})\times\frac{s}{d} \end{aligned} $$ trick9：当枚举的是所有因数的时候，枚举 $\dfrac{n}{d}$ 和 $d$ 是一样的。 $$ ans=\sum\limits_{s=1}^{\min(n,m)}s\times sum(\frac{n}{s})\times sum(\frac{m}{s})\sum\limits_{d=1,d\mid s}\mu(d)\times d $$ 后面的求和似乎是个积性函数，可以筛出来然后求前缀和。然后同样整除分块。
Code
写的是两个整除分块的做法。
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
using std::cin;
using std::cout;
using std::endl;
const int maxn = 1e7 + 10, mod = 20101009, inv_2 = 10050505;
typedef long long ll;
inline void upd(ll &x, ll y) { x += y, x %= mod; }
int pri[maxn], mu[maxn], tot, vis[maxn];
ll f[maxn];
void init(const int n = 1e7)
{
mu[1] = 1;
for (int i = 2; i <= n; i++)
{
if (!vis[i])
pri[++tot] = i, mu[i] = -1;
for (int j = 1; j <= tot && (ll)i * pri[j] <= n; j++)
{
vis[pri[j] * i] = 1;
if (i % pri[j] == 0)
{
mu[i * pri[j]] = 0;
break;
}
mu[i * pri[j]] = -mu[i];
}
}
for (int i = 1; i <= n; i++)
f[i] = f[i - 1], upd(f[i], 1ll * i * i % mod * (mu[i] + mod) % mod);
}
ll S(ll x) { return x * (x + 1) / 2 % mod; }
ll solve(int n, int m)
{
ll res = 0;
for (int i = 1, j; i <= std::min(n, m); i = j + 1)
{
j = std::min(n / (n / i), m / (m / i));
upd(res, (f[j] - f[i - 1] + mod) % mod * S(n / i) % mod * S(m / i) % mod);
}
return res;
}
int main()
{
int n, m;
cin >> n >> m;
ll ans = 0;
init();
for (int i = 1, j; i <= std::min(n, m); i = j + 1)
{
j = std::min(n / (n / i), m / (m / i));
upd(ans, solve(n / i, m / i) * (i + j) % mod * inv_2 % mod * (j - i + 1) % mod);
}
cout << ans << endl;
}