【bzoj2956】模积和

　求∑∑((n mod i)*(m mod j))其中1<=i<=n,1<=j<=m,i≠j。

若没有i不等于j，那么只需求sigma(n%i) * sigma(m%i)。这个可以发现商的种类只有sqrt(n)级别。枚举商发现余数是等差数列。

对于i = j，和刚才类似，考虑在n / i与m / i不变时，由两个等差数列。(s1 + i * p) * (s2 + i * q) 拆开即可。

#include <iostream>
#include <algorithm>
#include <cstdio>
#include <string>
#include <cstring>
#include <ctime>
#define Rep(i, x, y) for (int i = x; i <= y; i ++)
#define RepE(i, x) for (int i = pos[x]; i; i = g[i].nex)
using namespace std;
typedef long long LL;
const int mod = 19940417;
int n, m;
LL ans;
LL Calc(int n) {
    LL ret = 0;
    for (int i = 1; i <= n; ) {
        int t = n / i, lt = n / t, st = n % lt, ed = n % i;
        (ret += LL(ed + st) * ((ed - st) / t + 1) / 2) %= mod;
        i = lt + 1;
    }
    return ret;
}
LL G(LL x) {
    LL o = x * (x + 1) / 2;
    if (o % 3 == 0) o = ((o / 3) % mod) * (2*x+1);
    else o = o % mod * (2*x + 1) / 3;
    return o % mod;
}
int main()
{
    scanf ("%d%d", &n, &m);
    if (n < m) swap(n, m);
    ans = Calc(n) * Calc(m) % mod;
    for (int i = 1, lt; i <= m; i = lt + 1) {
        LL ta = n / i, tb = m / i;
        if (i <= m) lt = min(n/ta, m/tb);
        else lt = n / ta;
        LL s1 = n % lt, e1 = n % i, s2 = m % lt, e2 = m % i, l = (e1 - s1) / ta + 1;
        LL ret = ((l * s1 % mod) * s2 + (l*(l-1)/2 % mod) * ((ta*s2 + tb*s1) % mod) + (G(l-1) * ta % mod) * tb) % mod;
        ans -= ret;
    }
    ans = ((ans % mod) + mod) % mod;
    printf ("%lld\n", ans);
    return 0;
}