题意
求 \(\sum\limits_{i = 1}^N \sum\limits_{j = 1}^M lcm (i, j)\)
Solution
易知，原式
\[\sum\limits_{i = 1}^N \sum\limits_{j = 1}^M \frac{ij}{\gcd (i, j)}\] 枚举 \(\gcd (i, j)\) ，且将 \(d\) 提出来得\[\sum\limits_{d = 1}^{\min (N, M)} d \sum\limits_{i = 1}^{\left\lfloor\frac{N}{d}\right\rfloor} \sum\limits_{j = 1}^{\left\lfloor\frac{M}{d}\right\rfloor} ij[(i, j) = 1]\] 将公式 \(\sum\limits_{k | n} \mu(k) = [n = 1]\) 代入，得\[\sum\limits_{d = 1}^{\min (N, M)} d \sum\limits_{i = 1}^{\left\lfloor\frac{N}{d}\right\rfloor} \sum\limits_{j = 1}^{\left\lfloor\frac{M}{d}\right\rfloor} ij \sum\limits_{k | (i, j)} \mu(k)\] 套路枚举 \(k\) ，得\[\sum\limits_{d = 1}^{\min (N, M)} d \sum\limits_{k = 1}^{\min (\left\lfloor\frac{N}{d}\right\rfloor, \left\lfloor\frac{M}{d}\right\rfloor)} \mu(k) \sum\limits_{i = 1}^{\left\lfloor\frac{N}{d}\right\rfloor} \sum\limits_{j = 1}^{\left\lfloor\frac{M}{d}\right\rfloor} ij [k | (i, j)]\] 那么 \(ij\) 存在贡献时其必定是 \(k\) 的倍数，故\[\sum\limits_{d = 1}^{\min (N, M)} d \sum\limits_{k = 1}^{\min (\left\lfloor\frac{N}{d}\right\rfloor, \left\lfloor\frac{M}{d}\right\rfloor)} \mu(k) \sum\limits_{ki = 1}^{\left\lfloor\frac{N}{d}\right\rfloor} \sum\limits_{kj = 1}^{\left\lfloor\frac{M}{d}\right\rfloor} k^2 ij\] 将 \(k\) 提出，得\[\sum\limits_{d = 1}^{\min (N, M)} d \sum\limits_{k = 1}^{\min (\left\lfloor\frac{N}{d}\right\rfloor, \left\lfloor\frac{M}{d}\right\rfloor)} k^2 \mu(k) ( \sum\limits_{i = 1}^{\left\lfloor\frac{N}{kd}\right\rfloor} i) (\sum\limits_{j = 1}^{\left\lfloor\frac{M}{kd}\right\rfloor} j)\] 那么就可以预处理 \(\sum\limits_{k = 1}^n k^2 \mu(k)\) ，后面的用整除分块就好了Code
#include#include #include #define MOD 20101009using namespace std;typedef long long LL;const int MAXN = 1e07 + 10;int prime[MAXN];int vis[MAXN]= {0};int pcnt = 0;int mu[MAXN]= {0};LL sum[MAXN]= {0};const int MAX = 1e07;void prime_Acqu () { mu[1] = 1; for (int i = 2; i <= MAX; i ++) { if (! vis[i]) { prime[++ pcnt] = i; mu[i] = - 1; } for (int j = 1; j <= pcnt && i * prime[j] <= MAX; j ++) { vis[i * prime[j]] = 1; if (! (i % prime[j])) break; mu[i * prime[j]] = - mu[i]; } } for (int i = 1; i <= MAX; i ++) sum[i] = (sum[i - 1] + 1ll * i * 1ll * i % MOD * mu[i] % MOD) % MOD;}int N, M;inline LL calc (int n) { LL fn = (LL) n; return (fn * (fn + 1) >> 1) % MOD;}LL Solve () { LL ans = 0; int limit = min (N, M); for (int d = 1; d <= limit; d ++) { LL total = 0; int minlim = min (N / d, M / d); for (int l = 1, r; l <= minlim; l = r + 1) { r = min ((N / d) / ((N / d) / l), (M / d) / ((M / d) / l)); total = (total + (sum[r] - sum[l - 1] + MOD) % MOD * calc (N / d / l) % MOD * calc (M / d / l) % MOD) % MOD; } ans = (ans + (LL) (d) * total % MOD) % MOD; } return ans;}int main () { prime_Acqu (); scanf ("%d%d", & N, & M); LL ans = Solve (); cout << ans << endl; return 0;}/*4 5*/