P6825 「EZEC-4」求和

P6825 「EZEC-4」求和

求：
$$\sum _{i=1}^n\sum _{j=1}^n\gcd (i,j{)}^{i+j}$$
先化简原式，有：
$$\begin{array}{rlr}& \sum _{i=1}^n\sum _{j=1}^n\gcd (i,j{)}^{i+j}& \\ & =\sum _{d=1}^n\sum _{i=1}^n\sum _{j=1}^n{d}^{i+j}[\gcd (i,j)=d]& \\ & =\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }{d}^{d\cdot (i+j)}[\gcd (i,j)=1]& \\ & =\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }{d}^{d\cdot (i+j)}\sum _{p|\gcd (i,j)}\mu (p)& \\ & =\sum _{d=1}^n\sum _{p=1}^{\lfloor \frac{n}{d}\rfloor }\mu (p)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }{d}^{d\cdot (i+j)}[p|i][p|j]& \\ & =\sum _{d=1}^n\sum _{p=1}^{\lfloor \frac{n}{d}\rfloor }\mu (p)\sum _{i=1}^{\lfloor \frac{n}{pd}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{pd}\rfloor }{d}^{dp\cdot (i+j)}& \\ & =\sum _{d=1}^n\sum _{p=1}^{\lfloor \frac{n}{d}\rfloor }\mu (p)\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{T}\rfloor }{d}^{T\cdot (i+j)}& T=dp\\ & =\sum _{d=1}^n\sum _{p=1}^{\lfloor \frac{n}{d}\rfloor }\mu (p)\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }(d^T{)}^i\sum _{j=1}^{\lfloor \frac{n}{T}\rfloor }(d^T{)}^j& \\ & =\sum _{d=1}^n\sum _{p=1}^{\lfloor \frac{n}{d}\rfloor }\mu (p){\left(\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }(d^T{)}^i\right)}^2& \end{array}$$

再看后面那部分，令 $x={\frac{T}{p}}^T$。

设 $S_n$ 为等比数列前 $n$ 项和，则有：
$$S_n=S_{\lfloor \frac{n}{2}\rfloor }(1+a^{\lfloor \frac{n}{2}\rfloor })+a^n[2\nmid n]$$
即可 $\mathcal{O}(\log n)$ 求解。

时间复杂度 $\mathcal{O}(n\log n)$。

#include <bits/stdc++.h>

using namespace std;

#define ll long long

#define pll pair<ll, ll>

#define mpr make_pair

#define fi first

#define se second

const int _ = 1.5e6 + 10;

int mod;

int read()
{
    int x = 0, f = 1;
    char c = getchar();
    while (c < '0' || c > '9')
    {
        if (c == '-')
            f = -1;
        c = getchar();
    }
    while (c >= '0' && c <= '9')
    {
        x = x * 10 + c - '0';
        c = getchar();
    }
    return x * f;
}

int qpow(int x, int y)
{
    int res = 1;
    while (y)
    {
        if (y & 1)
            res = 1ll * res * x % mod;
        x = 1ll * x * x % mod;
        y >>= 1;
    }
    return res;
}

pll qsum(int a, int n)
{
    if (n == 0)
        return mpr(1, 1);
    if (n == 1)
        return mpr(a, a);
    pll tmp = qsum(a, n / 2);
    if (n & 1)
        return mpr((tmp.fi * (tmp.se + 1) % mod + tmp.se * tmp.se % mod * a % mod) % mod, tmp.se * tmp.se % mod * a % mod);
    else
        return mpr(tmp.fi * (tmp.se + 1) % mod, tmp.se * tmp.se % mod);
}

int T, n, m, k;

bitset<_> vis;

int cnt, pri[_], mu[_];

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++)
    {
        if (!vis[i])
            pri[++cnt] = i, mu[i] = -1;
        for (int j = 1; j <= cnt && i * pri[j] <= n; j++)
        {
            vis[i * pri[j]] = 1;
            if (i % pri[j] == 0)
                break;
            mu[i * pri[j]] = -mu[i];
        }
    }
}

int solve(int n)
{
    int ans = 0;
    for (int d = 1; d <= n; d++)
    {
        int a = qpow(d, d);
        int b = a;
        for (int p = 1; p <= n / d; p++)
        {
            pll tmp = qsum(b, n / (d * p));
            ans = (ans + tmp.fi * tmp.fi % mod * (mu[p] + mod) % mod) % mod;
            b = 1ll * b * a % mod;
        }
    }
    return ans;
}

signed main()
{
    T = read();
    init(_ - 10);
    while (T--)
    {
        n = read(), mod = read();
        printf("%d\n", solve(n));
    }
    return 0;
}