Prime prefix statistics
Find the sum of all prime numbers within \ (n \) to the power of \ (c \).
Consider the Ehrlich sieve, enumerate all multiples of it with a prime number at a time, and sieve out all illegal numbers.
Some definitions and explanations
- \(P {\ min} (x) \): the minimum prime factor of \ (x \), \ (P {\ min} (1) = + \ infty \), if \ (x \) is not a prime number, then \ (P {\ min} (x) \ Le \ sqrt {x} \).
- \The number of prime numbers within (m: \ left\lfloor \sqrt{n}\right\rfloor \).
- \(p_k \): the \ (K \) prime number.
- \(P_k \): the set composed of the first \ (K \) prime numbers.
- \(S(n, k) \): the number within \ (n \) is the set of the remaining numbers after deleting the illegal number with the first \ (K \) prime number, that is \ (S(n, k) = \{x \in \mathbb{N}^* | x \le n \wedge p_{\min}(x) \ge p_k\}\). At the same time, we note that for \ (S(n, m) \) is all prime numbers that have not been traversed.
- \(D(n, k) \): (D(n, k) = S(n, k - 1) - S(n, k) \) the set of numbers deleted by the prime number \ (k \) within \ (n \).
- \(f(n, k) \): \ (\ displaystyle \ sum {x \ in S (n, K)} x ^ C \), so the answer to the question is \ (f(n, m) - 1 + \displaystyle \sum_{x \in P_m} x^c\).
Derivation process
\(S(n, k) = S(n, k - 1) - p_k S(\left\lfloor \frac{n}{p_k}\right\rfloor , k - 1)\)
That is, consider the number deleted by \ (p_k \), which can be proved by the fact that \ (p_k s (\ left \ lfloor \ frac{n}{p_k} \ right \ rfloor, K - 1) \) and \ (D(n, k) \) are bijective.
Therefore, we can get the formula of \ (f \):
\(f(n, k) = f(n, k - 1) - p_k^c f(\left\lfloor\frac{n}{p_k}\right\rfloor, k - 1)\)
If we directly use the above formula to calculate each prime number, there will be \ (\ frac {\ sqrt{n} {\ log n} \) prime numbers in total. According to the division and blocking, we actually need to calculate \ (\ sqrt{n} \), so the complexity is \ (\ mathcal O(\frac{n}{\log n}) \)!
We can consider optimization. We will notice that for \ (k > m \), \ (S(n, k) = S(n, k - 1) - \{p_k \} \), because \ (S \) is still changing, we have to calculate the values of those \ (k > m \). If we let \ (S \) not change, will there be some optimization?
We can define something new:
- \(T(n, k) = S(n, k) +P_k\)
- \(g(n, k) = f(n, k) + \displaystyle \sum_{x \in P_m} x^c \), and the final answer is \ (g(n, m) - 1 \)
So we can use the formula of \ (S \) to get:
\(T(n, k) = T(n, k - 1) - p_k \left(T(\left\lfloor \frac{n}{p_k}\right\rfloor , k - 1) - P_{k - 1}\right) + \{p_k\}\)
We also note that \ (T(p_k - 1, k - 1) \) is \ (P_{k - 1} + {1} \), so we can also get the formula of \ (g \):
\(g(n, k) = g(n, k - 1) - \left(g(\left\lfloor \frac{n}{p_k}\right\rfloor , k - 1) - g(p_k - 1, k - 1)\right)\)
As for \ (g(n, 0) = \sum_{i \le n} i^c \), which can be obtained by Lagrange interpolation in \ (\ mathcal O(c) \).
Then the complexity is \ (\ mathcal o (\ dfrac {n ^ {3 / 4} {\ log n}) \), which is not easy to prove.
Template tea: Luogu5493 prime prefix statistics
#define in read<int>()
#define lin read<ll>()
#define rep(i, x, y) for(int i = (x); i <= (y); i++)
#define per(i, x, y) for(int i = (x); i >= (y); i--)
using namespace std;
using ll = long long;
template < typename T > void chkmax(T &x, const T &y) { x = x > y ? x : y; }
template < typename T > void chkmin(T &x, const T &y) { x = x < y ? x : y; }
const int N = 1e6 + 10;
const int inf = 0x7fffffff;
int mod, C;
ll lim, n;
ll g0[N], g1[N];
// g0[i] is g(i), g1[i] is g(n / i)
ll qp(ll x, int t) { x %= mod; ll res = 1; for(; t; t = t >> 1, x = x * x % mod) if(t & 1) res = res * x % mod; return res; }
namespace calcg {
const int N = 20;
int y[N], fac[N], ifac[N], pre[N], suf[N], m;
void init() {
m = C + 2;
fac[0] = 1; rep(i, 1, m) fac[i] = 1ll * fac[i - 1] * i % mod;
ifac[m] = qp(fac[m], mod - 2);
per(i, m - 1, 0) ifac[i] = 1ll * ifac[i + 1] * (i + 1) % mod;
rep(i, 1, m) y[i] = qp(i, C), y[i] = (y[i] + y[i - 1]) % mod;
}
int calc(ll v) {
v %= mod;
pre[0] = 1; rep(i, 1, m) pre[i] = 1ll * pre[i - 1] * (v - i + mod) % mod;
suf[m + 1] = 1; per(i, m, 1) suf[i] = 1ll * suf[i + 1] * (v - i + mod) % mod;
ll res = 0;
rep(i, 1, m) {
res += 1ll * pre[i - 1] * suf[i + 1] % mod * ifac[i - 1] % mod
* ifac[m - i] % mod * y[i] % mod * (m - i & 1 ? mod - 1 : 1) % mod;
res %= mod;
} return res;
}
}
using calcg :: calc;
int main() {
n = lin, C = in, mod = in, lim = sqrtl(n);
calcg :: init();
rep(i, 1, lim) g0[i] = (calc(i) - 1) % mod, g1[i] = (calc(n / i) - 1) % mod;
rep(i, 2, lim) {
g0[i] = (g0[i] + mod) % mod; if(g0[i] == g0[i - 1]) continue;
ll pc = qp(i, C), pg = g0[i - 1];
int l1 = lim / i, l2 = min(n / i / i, lim); chkmin(l1, l2);
rep(j, 1, l1) g1[j] = (g1[j] - pc * (g1[j * i] - pg)) % mod;
ll tt = n / i;
if(tt < inf) { // Kachang
int t = tt;
rep(j, l1 + 1, l2) g1[j] = (g1[j] - pc * (g0[t / j] - pg)) % mod;
} else {
ll t = tt;
rep(j, l1 + 1, l2) g1[j] = (g1[j] - pc * (g0[t / j] - pg)) % mod;
}
ll st = 1ll * i * i;
per(j, lim, st) g0[j] = (g0[j] - pc * (g0[j / i] - pg)) % mod;
}
ll ans = 0;
rep(i, 1, lim) {
g1[i] = (g1[i] + mod) % mod;
ans += 1ll * i * i % mod * g1[i] % mod; ans %= mod;
} printf("%lld\n", ans);
return 0;
}
Of course, this algorithm can make a wave of magic changes, and then you can do it LOJ6028 "from CommonAnts" prime count II Well, code.
Min_25 sieve
Give an integral function \ (F \):
- \(F(p) \) is a simple polynomial about \ (p \).
- \(F(p^k) \) can be obtained quickly.
- Find \ (\ sum {I = 1} ^ n f (I) \).
Consider the number of \ (P {\ min} = k \) in the order of \ (K \) from large to small, then record \ (S(n, k) = \{x \in \mathbb{N}^* | x \le n \wedge p_{\min}(x) \ge p_k \} \), remember \ (f(n, k) = \sum_{x \in S(n, k)} F(x) \), and the final answer is \ (f(n, 1) \).
So there is \ (\ displaystyle f(n, k) = \sum_{i = k}^{p_i \le n} \sum_{v} F(p_i ^ v)f\left(\left\lfloor\frac{n}{p_i ^ v}\right\rfloor, k+ 1\right) \), when \ (p_i \ge \sqrt{n} \), the value of the latter formula is \ (F(p_i) \), so it can be solved with the prime prefix and prefix preprocessed before.
The complexity is \ (\ mathcal o (\ dfrac {n ^ {3 / 4} {\ log n}) \), which is unlikely to prove.
Better than the above \ (\ downarrow \)
#include <bits/stdc++.h>
#define in read<int>()
#define lin read<ll>()
#define rep(i, x, y) for(int i = (x); i <= (y); i++)
#define per(i, x, y) for(int i = (x); i >= (y); i--)
using namespace std;
using ll = long long;
const int N = 1e6 + 10;
const int mod = 1e9 + 7;
const int inv2 = (mod + 1) >> 1;
const int inv3 = (mod + 1) / 3;
ll n, lim, sp1[N], sp2[N], g1[N], g2[N], w[N];
int ind1[N], ind2[N], tot;
int pri[N / 10], pnum;
bool v[N];
void shai(int l) {
rep(i, 2, l) {
if(!v[i]) {
pri[++pnum] = i;
sp1[pnum] = (sp1[pnum - 1] + i) % mod;
sp2[pnum] = (sp2[pnum - 1] + 1ll * i * i) % mod;
}
rep(j, 1, pnum) {
if(pri[j] * i > l) break;
v[pri[j] * i] = true; if(i % pri[j] == 0) break;
}
}
}
int S(ll x, int y) {
if(x < pri[y]) return 0;
ll res = 0;
res += g2[x <= lim ? ind1[x] : ind2[n / x]] - sp2[y - 1];
res -= g1[x <= lim ? ind1[x] : ind2[n / x]] - sp1[y - 1];
res = (res % mod + mod) % mod;
rep(k, y, pnum) {
if(1ll * pri[k] * pri[k] > x) break;
ll cur = pri[k];
for(int c = 1; cur <= x; cur *= pri[k], c++) {
res += cur % mod * ((cur - 1) % mod) % mod * (S(x / cur, k + 1) + (c != 1)) % mod;
res %= mod;
} // Enumerate the number of calculations, when c= When 1, it should include 1. When c = 1, we have calculated it.
} return res % mod;
}
int main() {
n = lin, lim = sqrtl(n);
shai(lim * 2); // Sift a little more
for(ll i = 1; i <= n; i++) {
++tot, w[tot] = (n / i); ll v = w[tot] % mod;
g1[tot] = (v * (v + 1) / 2 - 1) % mod;
g2[tot] = (v * (v + 1) % mod * (v * 2 + 1) % mod * inv2 % mod * inv3 % mod - 1) % mod;
if(w[tot] <= lim) ind1[w[tot]] = tot; else ind2[n / w[tot]] = tot;
i = n / (n / i);
}
rep(i, 1, pnum) {
int p = pri[i]; ll x1 = sp1[i - 1], x2 = sp2[i - 1];
rep(x, 1, tot) {
if(1ll * p * p > w[x]) break;
int y = w[x] / p <= lim ? ind1[w[x] / p] : ind2[n / (w[x] / p)];
g1[x] = (g1[x] - 1ll * p * (g1[y] - x1) % mod + mod) % mod;
g2[x] = (g2[x] - 1ll * p * p % mod * (g2[y] - x2) % mod + mod) % mod;
}
}
printf("%d\n", (S(n, 1) + 1) % mod);
return 0;
}