class="markdown_views prism-atom-one-light"> Description
求∑∑((n mod i)*(m mod j))其中1<=i<=n,1<=j<=m,i≠j。

---

Sample Input
3 4

---

Sample Output
1

---

这道题数论分块，我才接触。。。
设min(n, m)为o把那个式子划一下变成
∑(n - [n / i] * i) (1<=i<=n) * ∑(m - [m / i] * i) (1<=i<=m) - ∑(n - [n / i] * i) * (m - [m / i] * i) (1<=i<=o)
(n * n - ∑([n / i] * i) (1<=i<=n)) * (m * m - ∑([m / i] * i)(1<=i<=m)) - n * m * o - ∑(n * [m / i] * i) (1<=i<=o) - ∑(m * [n / i] * i) (1<=i<=o) + ∑([n / i] * i * [m / i] * i) (1<=i<=o)
设∑([n / i] * i) (1<=i<=n)为s1，∑([m / i] * i)(1<=i<=m)为s2。
设∑(n * [m / i] * i) (1<=i<=o)为s3，∑(m * [n / i] * i) (1<=i<=o)为s4。
设∑([n / i] * [m / i] * i * i) (1<=i<=o)为s5。
i * i的前缀和这部分可以靠(i * (i + 1) * (2 * i + 1)) / 6 线性求，于是求个逆元。
不过我离线搞了，但实际上也可以不求逆元，只是我嫌麻烦。。。
s1，s2，s3，s4，s5可以数论分块搞，其他都是O(1)的。

---

#include #include using namespace std; typedef long long LL; LL _min(LL x, LL y) {return x < y ? x : y;} LL _max(LL x, LL y) {return x > y ? x : y;} const LL mod = 19940417; const LL phi6 = 3323403; LL s1, s2; int main() { LL n, m; scanf("%lld%lld", &n, &m); s1 = s2 = 0; for(LL l = 1, r; l <= n; l = r + 1) { r = n / (n / l); LL pp = (r + l) * (r - l + 1) / 2; pp %= mod; (s1 += (pp * (n / l)) % mod) %= mod; } for(LL l = 1, r; l <= m; l = r + 1) { r = m / (m / l); LL pp = (r + l) * (r - l + 1) / 2; pp %= mod; (s2 += (pp * (m / l)) % mod) %= mod; } LL nm = (n * m) % mod, uu = (nm * nm) % mod; (uu -= ((n * n) % mod * s2) % mod) % mod; (uu -= ((m * m) % mod * s1) % mod) % mod; (uu += (s1 * s2) % mod) %= mod; LL oo = _min(n, m); s1 = s2 = 0; for(LL l = 1, r; l <= oo; l = r + 1) { r = _min(n / (n / l), oo); LL pp = (r + l) * (r - l + 1) / 2; pp %= mod; (s1 += (pp * (n / l)) % mod) %= mod; } for(LL l = 1, r; l <= oo; l = r + 1) { r = _min(m / (m / l), oo); LL pp = (r + l) * (r - l + 1) / 2; pp %= mod; (s2 += (pp * (m / l)) % mod) %= mod; } (uu -= (nm * oo) % mod) % mod; (uu += (n * s2) % mod); (uu += (m * s1) % mod); LL l1 = 1, r1 = 1, l2 = 1, r2 = 1, pp = 0; while(1) { LL ll = _max(l1, l2) - 1, rr = _min(r1, r2); LL ls = (((ll * 2 + 1) * (ll + 1) % mod * ll) % mod * phi6) % mod; LL rs = (((rr * 2 + 1) * (rr + 1) % mod * rr) % mod * phi6) % mod; (pp += ((n / l1) * (m / l2)) % mod * (rs - ls) % mod) %= mod; if(r1 == oo && r2 == oo) break; if(r1 == r2) { l1 = r1 + 1; r1 = _min(n / (n / l1), oo); l2 = r2 + 1; r2 = _min(m / (m / l2), oo); } else if(r1 < r2) { l1 = r1 + 1; r1 = _min(n / (n / l1), oo); } else { l2 = r2 + 1; r2 = _min(m / (m / l2), oo); } } (uu -= pp) %= mod; printf("%lld ", (uu + mod) % mod); return 0; }