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;
}

浙公网安备 33010602011771号