Min25 筛

我们将要通过 Min25 筛解决积性函数前缀和的问题。

\(P_k\) 为第 \(k\) 个质数,我们将所有数分成质数和合数进行讨论。并且我们不妨将 \(p^k(p^k - 1)\) 拆成 \(p^{2k} - p^k\),接下来只讨论对于 \(p^k\) 的求和。


\(S(n, j) = \sum_{i = 2}^nf(i)[\underset{p | i}{\min} p > P_j]\)

质数

\(q(n) = \sum_{i = 2}^n f(i)[i \in Prime]\),那么 \(S(n, j)\) 在质数处的和可以容斥为 \(q(n) - \sum_{i = 1}^j f(P_i)\)

合数

枚举一个数最小的质因子和它的次数,应为 \(\sum_{k>j}\sum_{e = 1}^{P_k^e \leq n}f(P_k^e)(S(\frac{n}{P_k^e}, k) + [e > 1])\)

综上,\(S(n, j) = q(n) - \sum_{i = 1}^j f(P_i) + \sum_{k>j}\sum_{e = 1}^{P_k^e \leq n}f(P_k^e)(S(\frac{n}{P_k^e}, k) + [e > 1])\)


现在考虑如何求 \(q(n)\)

设完全积性函数 \(f'\) 在质数处取值与 \(f\) 相同。显然 \(f'(p^k) = p^k\) 就是一种选择。

\(g(n, j) = \sum_{i = 2}^n f'(i)[i \in Prime \vee \underset{p | i}{\min} p > P_j]\),那么 \(q(n) = g(n, |P|)\)。考虑计算 \(g(n, j - 1) - g(n, j)\)

\(P_j^2 > n\),显然 \(g(n, j - 1) = g(n, j)\)

\(P_j^2 \leq n\),此时多出的贡献一定包含质因子 \(P_j\),因此为 \(g(n, j) = g(n, j - 1) - f'(P_j)(g(\frac{n}{P_j}, j - 1) - \sum_{i = 1}^{j - 1}f'(P_i))\)

会发现用到的下标类似 \(\lfloor \frac{n}{k} \rfloor\),取值只有 \(2\sqrt n\) 个。让 \(id_1\) 记录 \(k \leq \sqrt n\) 的下标,\(id_2\) 记录 \(k > \sqrt n\) 的下标。


#include <bits/stdc++.h>
using namespace std;
#define rep(i, a, b) for (int i = (a); i <= (b); i ++)
#define fro(i, a, b) for (int i = (a); i >= b; i --)
#define INF 0x3f3f3f3f
#define eps 1e-6
#define lowbit(x) (x & (-x))
#define reg register
#define IL inline
typedef long long LL;
typedef std::pair<int, int> PII;
#define int long long
inline int read() {
    int x = 0, f = 1;
    char ch = getchar();
    while (ch < '0' || ch > '9') { if (ch == '-') f = -1; ch = getchar(); }
    while (ch >= '0' && ch <= '9') { x = (x << 1) + (x << 3) + (ch ^ 48); ch = getchar(); }
    return x * f;
}

const int N = 200010, Mod = 1e9 + 7;

int n, m;
int w[N], id1[N], id2[N], g1[N], g2[N], sq;
int primes[N], cnt, st[N], s1[N], s2[N];

void prework() {
    int n = sq;
    for (int i = 2; i <= n; i ++) {
        if (!st[i]) primes[++ cnt] = i;
        for (int j = 1; primes[j] <= n / i; j ++) {
            st[i * primes[j]] = 1;
            if (i % primes[j] == 0) break;
        }
    }
    for (int i = 1; i <= cnt; i ++) {
        s1[i] = (s1[i - 1] + primes[i]) % Mod;
        s2[i] = (s2[i - 1] + primes[i] * primes[i] % Mod) % Mod;
    }
}

inline int getid(int x) {
    if (x <= sq) return id1[x];
    return id2[n / x];
}
inline int f1(int x) { return x % Mod * (x % Mod + 1) / 2 % Mod; }
inline int f2(int x) { return x % Mod * (x % Mod + 1) % Mod * (2 * x % Mod + 1) % Mod * 166666668 % Mod; }

int S(int x, int j) {
    if (primes[j] > x) return 0;
    int ans = (((g2[getid(x)] - g1[getid(x)]) % Mod - (s2[j] - s1[j]) % Mod) % Mod + Mod) % Mod;
    for (int k = j + 1; k <= cnt && primes[k] * primes[k] <= x; k ++) 
        for (int e = 1, sp = primes[k]; sp <= x; sp *= primes[k], e ++)
            ans = (ans + sp % Mod * (sp % Mod - 1) % Mod * (S(x / sp, k) + (e > 1)) % Mod) % Mod;
    return ans;
}

signed main() {
    n = read(), sq = sqrt((LL)n), prework();
    for (int l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l); w[++ m] = n / l;
        g1[m] = f1(w[m]) - 1, g2[m] = f2(w[m]) - 1;
        if (w[m] <= sq) id1[w[m]] = m;
        else id2[n / w[m]] = m;
    }
    for (int i = 1; i <= cnt; i ++) 
        for (int j = 1; j <= m && primes[i] * primes[i] <= w[j]; j ++) {
            g1[j] = (g1[j] - primes[i] * (g1[getid(w[j] / primes[i])] - s1[i - 1]) % Mod + Mod) % Mod;
            g2[j] = (g2[j] - primes[i] * primes[i] % Mod * (g2[getid(w[j] / primes[i])] - s2[i - 1]) % Mod + Mod) % Mod;
        }
    LL ans = (S(n, 0) + 1) % Mod;
    printf("%lld\n", ans);
    return 0;
}
posted @ 2026-03-17 09:07  はなこくん  阅读(6)  评论(0)    收藏  举报