Luogu P6825 「EZEC-4」求和

题目描述

https://www.luogu.com.cn/problem/P6825

简要题意：给定 $n,p$，求 $\sum_{i=1}^n\sum_{j=1}^n(i,j)^{(i,j)}$

$n\le 1.5\times 10^6$

Solution

经过简单的式子化简，我们可以得到 $ans=\sum_{T=1}^n\sum_{t|T}\mu(\frac{T}{t})f^2(t^T,\lfloor\frac{n}{T}\rfloor)$，其中 $f(x,y)=\sum_{i=1}^yx^i$

我们考虑直接枚举 $T,t$ 这个东西的时间复杂度是 $O(n\log n)$，$f(x,y)$ 显然是一个等比数列求和东西，我们直接套用等比数列求和的公式的话的时间复杂度为 $O(\log p)$，总的时间复杂度为 $O(n\log n\log p)$

但是我们注意到 $t^T$ 的前缀和最多才到 $\lfloor\frac{n}{T}\rfloor$ 项，如果我们计算这部分的时间复杂度是 $O(\log \lfloor\frac{n}{T}\rfloor)$ 应该能显著减少常数，同时我想起来等比数列前缀和(n项)有 $O(\log n)$ 的做法，然后我们把这个做法套上去，发现这样就过了

但其实 $O(\sum_{i=1}^n\sum_{j=1}^{\lfloor\frac{n}{i}\rfloor}\log(\lfloor\frac{n}{ij}\rfloor)=O(n\log n)$，感觉非常神奇

#include <iostream>
#define maxn 1500010
#define ll long long
using namespace std;

int p;
inline int add(int x, int y) { return (x += y) >= p ? x - p : x; }
inline int mul(int x, int y) { return 1ll * x * y % p; }
inline int add(initializer_list<int> lst) { int s = 0; for (auto t : lst) s = add(s, t); return s; }
inline int mul(initializer_list<int> lst) { int s = 1; for (auto t : lst) s = mul(s, t); return s; }
int pow_mod(int x, int n) {
    int s = 1;
    for (; n; n >>= 1, x = mul(x, x))
        if (n & 1) s = mul(s, x);
    return s;
}

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

int f(int a, int n) {
    if (n == 0) return 1;
    if (n & 1) return mul(f(a, n / 2), 1 + pow_mod(a, n / 2 + 1));
    else return add(mul(f(a, n / 2 - 1), 1 + pow_mod(a, n / 2 + 1)), pow_mod(a, n / 2));
}

int n;

void work() {
    cin >> n >> p; int ans = 0;
    for (int i = 1; i <= n; ++i) {
        int d = pow_mod(i, i), D = 1;
        for (int j = 1; j <= n / i; ++j) {
            D = mul(D, d); if (!mu[j]) continue;
            int _mu = mu[j] < 0 ? p - 1 : 1, t = f(D, n / (i * j)) - 1; 
            ans = add(ans, mul({ _mu, t, t }));
        }
    } cout << ans << "\n";
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr); cout.tie(nullptr);

    int T; cin >> T; init_isp(1500000);
    while (T--) work(); 
    return 0;
}