P3768:简单的数学题题解

解法

题意就是要求:

\[ ( \sum_{i = 1}^n \sum_{j = 1}^n ij gcd(i, j) ) \mod p \]

言简意赅的题面非常舒适。出现了\(gcd\)可以考虑往反演方面思考。先进行例行的套路:枚举\(gcd\)。在这里为了书写方便,我就把取模省去。

\[ \sum_{d = 1}^n d \sum_{i = 1}^n \sum_{j = 1}^n ij [gcd(i, j) = d] \]

然后在后面消掉\(d\)这个玩意,同时约去两个\(d\)就会产生\(d^2\),乘到外面即可:

\[ \sum_{d = 1}^n d^3 \sum_{i = 1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j = 1}^{\lfloor \frac{n}{d} \rfloor} ij [gcd(i, j) = 1] \]

发现后面的式子可以直接反演,考虑利用经典性质代换:

\[ \sum_{d = 1}^n d^3 \sum_{i = 1}^{\lfloor \frac{n}{d} \rfloor} i \sum_{j = 1}^{\lfloor \frac{n}{d} \rfloor} j \sum_{h = 1}^n \mu(h) [h|i][h|j] \]

调换式子位置:

\[ \sum_{d = 1}^n d^3 \sum_{h = 1}^n \mu(h) \sum_{i = 1}^{\lfloor \frac{n}{d} \rfloor} i[h|i] \sum_{j = 1}^{\lfloor \frac{n}{d} \rfloor} j[h|j] \]

除掉两个\(h\)出来:

\[ \sum_{d = 1}^n d^3 \sum_{h = 1}^n \mu(h) h^2 \sum_{i = 1}^{\lfloor \frac{n}{hd} \rfloor} i \sum_{j = 1}^{\lfloor \frac{n}{hd} \rfloor} j \]

后面两个好像可以用高中数列求和用\(O(1)\)的时间求出,设\(calc(n) = \sum_{i = 1}^n i = \frac{n(n + 1)}{2}\):

\[ \sum_{d = 1}^n d^3 \sum_{h = 1}^n \mu(h) h^2 calc(\lfloor \frac{n}{hd} \rfloor)^2 \]

设\(T = hd\)并化简:

\[ \sum_{d = 1}^n d^3 \sum_{T = 1}^n \mu(\lfloor \frac{T}{d} \rfloor) (\lfloor \frac{T}{d} \rfloor)^2 calc(\lfloor \frac{n}{T} \rfloor)^2 \\ \sum_{T = 1}^n T^2 calc(\lfloor \frac{n}{T} \rfloor)^2 \sum_{d|T} d \mu(\lfloor \frac{T}{d} \rfloor) \]

交换\(T\)和\(d\)位置时,要注意\(T = hd\)这个条件,故\(d|T\)需要被满足。我们发现后面这个式子:

\[ \sum_{d|T} d \mu(\lfloor \frac{T}{d} \rfloor) \]

这个式子是 Dirichlet 卷积的形式,简写为:

\[ id * \mu \]

对 Dirichlet 卷积熟悉的读者一眼可以看出来这两个函数的卷积为\(\varphi\)。这个结论一行就可以证完,如果不知道\( \epsilon = \mu * 1 \)怎么来的读者可以阅读这里,以下是证明:

\[ \varphi * 1 = id = id * \epsilon = id * \mu * 1 = id * \mu = \varphi \]

所以直接代换掉:

\[ \sum_{T = 1}^n calc(\lfloor \frac{n}{T} \rfloor)^2 \varphi(T) T^2 \]

已经很优秀了,用杜教筛搞下\(\varphi(T) T^2\)前缀和然后整除分块即可。

代码

// P3768.cpp
#include <bits/stdc++.h>
#define ll long long
using namespace std;

const ll MAX_N = 5e6 + 2000;
ll n, mod, phi[MAX_N], prime[MAX_N], tot, inv2, inv6;
bool vis[MAX_N];
unordered_map<ll, ll> ump;

ll getPreQuad(ll num)
{
    num %= mod;
    return num * (num + 1) % mod * inv2 % mod;
}

ll getPreThird(ll num)
{
    num %= mod;
    return num * (num + 1) % mod * (num + num + 1) % mod * inv6 % mod;
}

void phi_sieve()
{
    phi[1] = 1;
    for (ll i = 2; i < MAX_N; i++)
    {
        if (!vis[i])
            prime[++tot] = i, phi[i] = i - 1;
        for (ll j = 1; j <= tot && 1ll * prime[j] * i < MAX_N; j++)
        {
            vis[i * prime[j]] = true;
            if (i % prime[j] == 0)
            {
                phi[i * prime[j]] = phi[i] * prime[j] % mod;
                break;
            }
            phi[i * prime[j]] = phi[i] * phi[prime[j]] % mod;
        }
    }
    for (ll i = 2; i < MAX_N; i++)
        phi[i] = (phi[i] * i % mod * i % mod + phi[i - 1]) % mod;
}

ll quick_pow(ll bas, ll tim)
{
    ll ans = 1;
    while (tim)
    {
        if (tim & 1)
            ans = ans * bas % mod;
        bas = bas * bas % mod;
        tim >>= 1;
    }
    return ans;
}

ll sieve(ll num)
{
    if (num < MAX_N)
        return phi[num];
    if (ump[num])
        return ump[num];
    ll ans = getPreQuad(num);
    ans = ans * ans % mod;
    for (ll i = 2, gx; i <= num; i = gx + 1)
    {
        gx = num / (num / i);
        ll fact = (getPreThird(gx) - getPreThird(i - 1) + mod) % mod;
        ans = (ans - fact * sieve(num / i) % mod + mod) % mod;
    }
    return ump[num] = (ans + mod) % mod;
}

signed main()
{
    scanf("%lld%lld", &mod, &n);
    inv2 = quick_pow(2, mod - 2), inv6 = quick_pow(6, mod - 2);
    phi_sieve();
    ll ans = 0;
    for (ll i = 1, gx; i <= n; i = gx + 1)
    {
        gx = n / (n / i);
        ll fact = getPreQuad(n / i);
        fact = fact * fact % mod;
        ll bfact = (sieve(gx) - sieve(i - 1) + mod) % mod;
        ans = (ans + bfact * fact % mod) % mod;
    }
    printf("%lld", (ans + mod) % mod);
    return 0;
}

发布者

Kal0rona

江西师大附中全机房最弱

发表评论

电子邮件地址不会被公开。 必填项已用*标注