模积和 HYSBZ

题目
求(sum_{i = 1}^nsum_{j = 1}^m (n mod i) (m mod j)) 且(i≠j)

(sum_{i = 1}^nsum_{j = 1}^m(n - lfloorfrac{n}{i} floor * i)(m - lfloorfrac{m}{j} floor * j) - sum_{i = 1}^{min(n, m)}(n - lfloorfrac{n}{i} floor * i)(m - lfloorfrac{m}{i} floor * i))
求(f(n) = sum_{i = 1}^n(n - lfloorfrac{n}{i} floor * i) = n^2 - sum_{i = 1}^ni * lfloorfrac{n}{i} floor)
(sum_{i = 1}^{min(n, m)}(n - lfloorfrac{n}{i} floor * i)(m - lfloorfrac{m}{i} floor * i))
令(n ≤ m)
$sum_{i = 1}^{n}(n - lfloorfrac{n}{i} floor * i)(m - lfloorfrac{m}{i} floor * i) = sum_{i= 1} ^ n n * m - n * ilfloorfrac{m}{i} floor - m * i * lfloorfrac{n}{i} floor + i * i * lfloorfrac{n}{i} floor * lfloorfrac{m}{i} floor $

#include <iostream>
#include <cstdio>
#define ll long long
using namespace std;
const int mod = 19940417;
const int inv = 3323403; //6 在mod下的逆元
ll sum(ll l, ll r){
    return (r - l + 1) * (l + r) / 2 % mod;
}
ll sum2(ll x){
    return x * (x + 1) % mod * (2 * x + 1) % mod * inv % mod;
}
ll cal(ll n){
    ll ans = 0;
    for(ll l = 1, r; l <= n; l = r + 1){
        r = n / (n / l);
        ans = (ans + n * (r - l + 1) % mod - sum(l, r) * (n / l)) % mod;
    }
    return (ans + mod) % mod;
}
int main(){
    ll n, m;
    scanf("%lld%lld", &n, &m);
    ll ans = cal(n) * cal(m) % mod;
    if(n > m)swap(n, m);
    ll s1, s2, s3;
    for(ll l = 1, r; l <= n; l = r + 1){
        r = min(n / (n / l), m / (m / l));
        s1 = n * m % mod * (r - l + 1) % mod;
        s2 = (n / l) * (m / l) % mod * (sum2(r) - sum2(l - 1) + mod) % mod;
        s3 = ((n / l) * m % mod + (m / l) * n % mod) % mod * sum(l, r) % mod;
        ans = (ans - (s1 + s2 - s3) % mod + mod) % mod;
    }
    printf("%lld
", ans);
    return 0;
}