找素数
要求求出1e11以内的素数个数
这题有趣的是代码只求了g函数不用求S函数,试分析原因
g函数的意义是:前n个数中所有素数的贡献
Σ
i
=
1
n
[
i
∈
P
r
i
m
e
]
f
(
i
)
\Sigma_{i = 1}^{n}[i \in Prime]f(i)
Σi=1n[i∈Prime]f(i)
而S函数是:前n个数的贡献
Σ
i
=
1
n
f
(
i
)
\Sigma_{i = 1}^{n}f(i)
Σi=1nf(i)(即积性函数的前缀和)
提要求素数个数,可构造函数
f
(
x
)
=
[
x
∈
P
r
i
m
e
]
f(x) = [x \in Prime]
f(x)=[x∈Prime]则
S
(
n
,
0
)
S(n, 0)
S(n,0)即为所求答案。特别的,在合数部分
f
(
i
)
=
0
f(i) = 0
f(i)=0所以素数的贡献即所有数的贡献
S
(
n
,
0
)
=
g
(
n
,
0
)
S(n, 0) = g(n, 0)
S(n,0)=g(n,0)
#include <bits/stdc++.h>
using namespace std;
const int N = 1e6 + 10;
typedef long long ll;
ll n, sq;
ll prime[N], cnt, vis[N], sp[N];///sp(sum_prime)前p个质数的贡献和
ll w[N], tot, id1[N], id2[N], g[N];
ll fx(ll x) {
return !vis[x];//素数返回1计入贡献, 合数返回0
}
void sieve(ll n) {
vis[1] = 1;
for (int i = 2; i <= n; i++) {
if (!vis[i]) {
prime[++cnt] = i;
sp[cnt] = cnt;
}
for (int j = 1; j <= cnt && i * prime[j] <= n; j++) {
vis[i * prime[j]] = 1;
if (i % prime[j] == 0)break;
}
}
}
int main() {
cin >> n;
sq = sqrt(n);
sieve(sq);
for (ll l = 1, r; l <= n; l = r + 1) {
r = min(n, n / (n / l));
w[++tot] = n / l;
g[tot] = w[tot] - 1;
w[tot] = n / l;
if (w[tot] <= sq) id1[w[tot]] = tot;
else id2[n / w[tot]] = tot;
}
for (int j = 1; j <= cnt; j++) {
for (int i = 1; i <= tot && prime[j] * prime[j] <= w[i]; i++) {
ll tmp = w[i] / prime[j];
ll p = tmp <= sq ? id1[tmp] : id2[n / tmp];
g[i] = (g[i] - fx(prime[j]) * (g[p] - sp[j - 1]));
}
}
printf("%lld", g[1]);
}
类似的题
#188. 【UR #13】Sanrd
要求出
l
−
r
l-r
l−r范围内每个数的次大质因子之和,特别的,
f
(
p
i
)
=
0
,
p
i
f(p_i) = 0, p_i
f(pi)=0,pi表示质数。在该题中质数的贡献为0,则在公式
S
(
n
,
j
)
=
(
g
(
n
,
∣
P
∣
)
−
Σ
i
=
1
j
f
(
p
i
)
)
+
Σ
k
=
j
+
1
p
k
2
≤
n
Σ
e
=
1
p
k
e
≤
n
(
f
(
p
k
e
)
∗
(
S
(
⌊
n
p
k
e
⌋
,
k
)
+
[
e
>
1
]
)
S(n, j) = (g(n, |P|) - \Sigma_{i=1}^{j}f(p_i) )+ \Sigma_{k = j + 1}^{p_k^2 \leq n}\Sigma_{e=1}^{p_k^e\leq n}(f(p_k^e)*(S(\lfloor\frac{n}{p_k^e}\rfloor, k)+[e>1])
S(n,j)=(g(n,∣P∣)−Σi=1jf(pi))+Σk=j+1pk2≤nΣe=1pke≤n(f(pke)∗(S(⌊pken⌋,k)+[e>1])中前半部分为0.,所以我们不需要计算g
大佬们的代码都构造了一个
F
(
x
)
=
1
F(x) = 1
F(x)=1的函数, 先分析道理
loj 6053简单的函数
做题前的疑问:
1为什么所有素数(除了2)一定要用p - 1表示, 2用p+1表示。而不能用p^1表示。
2为什么两个函数要分开算
f
(
x
)
=
x
2
−
x
f(x) = x^2-x
f(x)=x2−x
g
(
x
)
=
x
2
g(x) = x^2
g(x)=x2
h
(
x
)
=
x
h(x) = x
h(x)=x
f
(
x
)
=
g
(
x
)
−
h
(
x
)
f(x) = g(x) - h(x)
f(x)=g(x)−h(x)三个函数中
f
(
x
)
f(x)
f(x)并不是积性函数, 把多项式拆成幂函数, 也就有了积性函数
#include <bits/stdc++.h>
using namespace std;
const int N = 1e6 + 10, mod = 1e9 + 7;
typedef long long ll;
ll n, sq;
ll vis[N], prime[N], cnt = 0, sp1[N], sp2[N];//g1(x) = x , g2(x) = 1
ll fastpow(ll base, ll pow) {
ll ans = 1;
while (pow) {
if (pow & 1)
ans = ans * base % mod;
base = base * base % mod;
pow >>= 1;
}
return ans;
}
ll inv2 = fastpow(2, mod - 2), inv6 = fastpow(6, mod - 2);
void sieve(ll n) {
vis[1] = 1;
for (int i = 2; i <= n; i++) {
if (!vis[i]) {
prime[++cnt] = i;
sp1[cnt] = (sp1[cnt - 1] + i) % mod;
sp2[cnt] = sp2[cnt - 1] + 1;
}
for (int j = 1; j <= cnt && i * prime[j] <= n; j++) {
vis[i * prime[j]] = 1;
if (i % prime[j] == 0)
break;
}
}
}
ll w[N], id1[N], id2[N], tot, g1[N], g2[N];
ll S(ll i, ll j) {
if (prime[j] > i)
return 0;
ll p = i <= sq ? id1[i] : id2[n / i];
ll ans = ((g1[p] - g2[p] + mod) % mod - (sp1[j] - sp2[j] + mod) % mod + mod) % mod;
if (i >= 2 && j == 0)
ans += 2;
for (int k = j + 1; k <= cnt && prime[k] * prime[k] <= i; k++) {
ll pe = prime[k];
for (int e = 1; pe <= i; e++, pe = pe * prime[k]) {
ll x = pe % mod;
ans = (ans + (prime[k] ^ e) % mod * (S(i / pe, k) + (e > 1)) % mod) % mod;
}
}
return ans;
}
inline ll cal_g1(ll x) {
return x;
}
inline ll cal_g2(ll x) {
return 1;
}
int main() {
cin >> n;
sq = sqrt(n);
sieve(sq);
//分块
for (ll l = 1, r; l <= n; l = r + 1) {
r = min(n / (n / l), n);
w[++tot] = n / l % mod;
//此处把所有数当成素数带入f(x), 因为题目条件只有素数和素数的次幂数能进行计算
//合数的计算并不方便, 为节省时间起见, 构造新函数F(x) = x ^ 2 - x
//因为合数最后会被筛去, 无论取值如何, 不影响结果
g1[tot] = (w[tot] * (w[tot] + 1) % mod * inv2 % mod - 1 + mod) % mod;
g2[tot] = (w[tot] - 1 + mod) % mod;
w[tot] = n / l;
if (w[tot] <= sq)
id1[w[tot]] = tot;
else
id2[n / w[tot]] = tot;
}
//计算g(n, cnt)
for (int j = 1; j <= cnt; j++) {
for (int i = 1; i <= tot && prime[j] * prime[j] <= w[i]; i++) {
ll temp = w[i] / prime[j];
ll p = temp <= sq ? id1[temp] : id2[n / temp];
g1[i] = (g1[i] - cal_g1(prime[j]) * (g1[p] - sp1[j - 1] + mod) % mod + mod) % mod;
g2[i] = (g2[i] - cal_g2(prime[j]) * (g2[p] - sp2[j - 1] + mod) % mod + mod) % mod;
}
}
printf("%lld", S(n, 0) + 1);
return 0;
}
#188. 【UR #13】Sanrd
要求出次大质因数
此题中比较纠结的点是g函数的适用范围,g函数表示1-n中素数的个数,在计算
g
(
n
,
j
)
=
g
(
n
,
j
−
1
)
−
F
(
p
j
)
∗
(
g
(
⌊
n
p
j
⌋
,
j
−
1
)
−
Σ
i
=
1
j
−
1
F
(
p
i
)
)
(
p
j
2
≤
n
)
g(n, j) = g(n, j - 1) - F(p_j) * (g(\lfloor\frac{n}{p_j}\rfloor, j - 1) - \Sigma_{i=1}^{j - 1}F(p_i)) (p_j^2 \leq n)
g(n,j)=g(n,j−1)−F(pj)∗(g(⌊pjn⌋,j−1)−Σi=1j−1F(pi))(pj2≤n)时,因为有条件
p
j
2
≤
n
p_j^2 \leq n
pj2≤n,可以保证
g
(
⌊
n
p
j
⌋
,
j
−
1
)
−
Σ
i
=
1
j
−
1
F
(
p
i
)
>
0
g(\lfloor\frac{n}{p_j}\rfloor, j - 1) - \Sigma_{i=1}^{j - 1}F(p_i) \gt 0
g(⌊pjn⌋,j−1)−Σi=1j−1F(pi)>0
使用时, 也要确保这个条件, 如此题计算S函数时, 设置
if (x < prime[i])continue;
的条件来保证
#include <bits/stdc++.h>
using namespace std;
#define Int register int
#define ll long long
#define MAXN 1000005
template <typename T> inline void read(T& t) { t = 0; char c = getchar(); int f = 1; while (c < '0' || c > '9') { if (c == '-') f = -f; c = getchar(); }while (c >= '0' && c <= '9') { t = (t << 3) + (t << 1) + c - '0'; c = getchar(); } t *= f; }
template <typename T, typename ... Args> inline void read(T& t, Args&... args) { read(t); read(args...); }
template <typename T> inline void write(T x) { if (x < 0) { x = -x; putchar('-'); }if (x > 9) write(x / 10); putchar(x % 10 + '0'); }
ll n, s, w[MAXN], g[MAXN];
int cnt, tot, id1[MAXN], id2[MAXN], vis[MAXN], prime[MAXN];
int& ID(ll x) { return x <= s ? id1[x] : id2[n / x]; }
void Euler(int up) {
for (Int i = 2; i <= up; ++i) {
if (!vis[i]) prime[++tot] = i;
for (Int j = 1; j <= tot && i * prime[j] <= up; ++j) {
vis[i * prime[j]] = 1;
if (i % prime[j] == 0) break;
}
}
}
ll S(ll x, ll y) {
if (prime[y] > x) return 0;
ll res = 0;
for (Int i = y + 1; i <= tot && 1ll * prime[i] * prime[i] <= x; ++i) {
for (ll pr = prime[i]; 1ll * pr * prime[i] <= x; pr *= prime[i]) {
int k = ID(x / pr);
if (x < prime[i])continue;
res += S(x / pr, i) + 1ll * prime[i] * (g[k] - i + 1);
//x会在不同的递归中被除以他的质因数, 若在某一次除法之后只剩一个质数
//则除数即为x的次大质因数,g[k] - i + 1表示大于等于p_i的素数个数(保证每个数只出现一次)
//g[k]表示1-w[k]中素数的个数
//g函数的隐式适用范围:
}
}
//printf("%lld %lld %lld\n", x, y, res);
return res;
}
ll Solve(ll N) {
cnt = tot = 0;
n = N, s = sqrt(n), Euler(s);
for (ll l = 1, r, v; l <= n; l = r + 1) r = (n / (v = n / l)), w[++cnt] = v, g[ID(v) = cnt] = v - 1;
for (Int i = 1; i <= tot; ++i)
for (Int j = 1; j <= cnt && 1ll * prime[i] * prime[i] <= w[j]; ++j) {
int k = ID(w[j] / prime[i]);
g[j] -= g[k] - (i - 1);
}
return S(n, 0);
}
signed main() {
ll l, r; read(l, r);
write(Solve(r) - Solve(l - 1)), putchar('\n');
return 0;
}