DDOSvoid's Blog

Luogu P5325 【模板】Min_25筛

发表于 分类于
本文字数： 1.9k 阅读时长 ≈ 2 分钟

题目描述

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

简要题意：定义积性函数 $f(x)$，且 $f(p^x)=(p^x)^2-p^x$，求 $\sum_{i=1}^nf(i)$

$n\le 10^{10}$

Solution

甚至不用化式子，直接使用 $Min\underline{}25$ 筛即可

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
#include <iostream>
#include <cmath>
#include <unordered_map>
#define ll long long
#define maxn 200010
using namespace std;

const int p = 1000000007;

ll pow_mod(ll x, ll n) {
    ll s = 1;
    for (; n; n >>= 1, x = x * x % p)
        if (n & 1) s = s * x % p;
    return s; 
}

ll n; 

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

ll id1[maxn], id2[maxn];
int get_id(ll x) { return x <= n / x ? id1[x] : id2[n / x]; }

ll g1[maxn], g2[maxn], w[maxn], inv2, inv6; int num;
void init_min25(ll n) {
    inv2 = pow_mod(2, p - 2); inv6 = pow_mod(6, p - 2);
    for (ll l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l); w[++num] = n / l; ll t = n / l % p;
        g1[num] = t * (t + 1) % p * inv2 % p - 1;
        g2[num] = t * (t + 1) % p * (2 * t + 1) % p * inv6 % p - 1; 
        if (w[num] <= n / w[num]) id1[w[num]] = num;
        else id2[n / w[num]] = num;
    }
    for (int i = 1; i <= cnt; ++i)
        for (int j = 1; j <= num && 1ll * pri[i] * pri[i] <= w[j]; ++j) {
            g1[j] = (g1[j] - pri[i] * (g1[get_id(w[j] / pri[i])] - sp1[i - 1])) % p;
            g2[j] = (g2[j] - 1ll * pri[i] * pri[i] % p * (g2[get_id(w[j] / pri[i])] - sp2[i - 1])) % p;
        }
}

ll S(ll n, int m) {
    if (pri[m] >= n) return 0;
    ll ans = (g2[get_id(n)] - g1[get_id(n)] - (sp2[m] - sp1[m])) % p;
    for (int i = m + 1; i <= cnt && 1ll * pri[i] * pri[i] <= n; ++i)
        for (ll x = 1, mul = pri[i]; mul <= n; mul *= pri[i], ++x)
            ans = (ans + mul % p * (mul % p - 1) % p * (S(n / mul, i) + (x > 1))) % p;
    return ans;
}

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

    cin >> n; init_isp(sqrt(n)); init_min25(n);
    cout << (S(n, 0) + 1 + p) % p << "\n";
    return 0;
}