P6810 「MCOI-02」Convex Hull 凸包

P6810 「MCOI-02」Convex Hull 凸包

推式子，不妨设 $n \le m$ ：
$$
\begin {aligned}
& \sum _ {i = 1} ^ n \sum _ {j = 1} ^ m d(i) d(j) d((i, j)) \\
= & \sum _ {x = 1} ^ {n} d(x) \sum _ {i = 1} ^ n \sum _ {j = 1} ^ m d(i) d(j) [(i, j) = x] \\
= & \sum _ {x = 1} ^ {n} d(x) \sum _ {i = 1} ^ {\left \lfloor \frac n x \right \rfloor} \sum _ {j = 1} ^ {\left \lfloor \frac m x \right \rfloor} d(ix) d(jx) [(i, j) = 1] \\
= & \sum _ {x = 1} ^ {n} d(x) \sum _ {i = 1} ^ {\left \lfloor \frac n x \right \rfloor} \sum _ {j = 1} ^ {\left \lfloor \frac m x \right \rfloor} d(ix) d(jx) \sum _ {y | i \wedge y | j} \mu(y) \\
= & \sum _ {x = 1} ^ {n} d(x) \sum _ {y = 1} ^ {\left \lfloor \frac n x \right \rfloor} \mu(y) \sum _ {i = 1} ^ {\left \lfloor \frac n {xy} \right \rfloor} \sum _ {j = 1} ^ {\left \lfloor \frac m {xy} \right \rfloor} d(ixy) d(jxy)
\end {aligned}
$$
令 $T = xy$ ，则：
$$
\begin {aligned}
& \sum _ {T = 1} ^ n \sum _ {x | T} d(x) \mu(\frac T x) \sum _ {i = 1} ^ {\left \lfloor \frac n T \right \rfloor} d(iT) \sum _ {j = 1} ^ {\left \lfloor \frac m T \right \rfloor} d(jT) \\
= & \sum _ {T = 1} ^ n d * \mu(T) \sum _ {i = 1} ^ {\left \lfloor \frac n T \right \rfloor} d(iT) \sum _ {j = 1} ^ {\left \lfloor \frac m T \right \rfloor} d(jT) \\
= & \sum _ {T = 1} ^ n \sum _ {i = 1} ^ {\left \lfloor \frac n T \right \rfloor} d(iT) \sum _ {j = 1} ^ {\left \lfloor \frac m T \right \rfloor} d(jT) \\
\end {aligned}
$$
注意到 $\displaystyle \sum _ {i = 1} ^ {\left \lfloor \frac n T \right \rfloor} d(iT)$ 与 $\displaystyle \sum _ {j = 1} ^ {\left \lfloor \frac m T \right \rfloor} d(jT)$ 是互相独立的，因此只有两层循环，复杂度为 $O(n \log n)$ 。

查看代码

#include <cstdio>
#include <algorithm>
using namespace std;
typedef long long LL;
const int N = 2e6 + 10;
int n, m, mod;
int cnt, p[N], d[N], low[N];
void init()
{
    d[1] = 1;
    low[1] = 1;
    for (int i = 2; i < N; i++)
    {
        if (!low[i])
        {
            p[++cnt] = i;
            d[i] = 2;
            low[i] = i;
        }
        for (int j = 1; j <= cnt && i * p[j] < N; j++)
        {
            if (i % p[j] == 0)
            {
                if (low[i] == i)
                    d[i * p[j]] = d[i] + 1;
                else
                    d[i * p[j]] = d[i / low[i]] * d[p[j] * low[i]];
                low[i * p[j]] = low[i] * p[j];
                break;
            }
            d[i * p[j]] = d[i] * d[p[j]];
            low[i * p[j]] = p[j];
        }
    }
}
int main()
{
    scanf("%d%d%d", &n, &m, &mod);
    if (n > m)
        swap(n, m);
    init();
    int res = 0;
    for (int i = 1; i <= n; i++)
    {
        int x = 0, y = 0;
        for (int j = 1; j <= n / i; j++)
            x = ((LL)x + d[i * j]) % mod;
        for (int j = 1; j <= m / i; j++)
            y = ((LL)y + d[i * j]) % mod;
        res = (res + (LL)x * y % mod) % mod;
    }
    printf("%d", res);
    return 0;
}