「Min_25 筛」学习笔记

Min_25 太神仙了!
据说这个筛法的复杂度是 Θ(n1ϵ)\Theta\left(n^{1 - \epsilon}\right) 的。(别问我 ϵ\epsilon 是啥,问就不知道
本质上其实是扩展埃筛。

感谢 chd 哥哥让我学会 Min_25 筛!(@Trisolaris

筛啥用的?

积性函数前缀和,不然用着干啥(
要求被筛的函数 ff 在质数 pp 处的值 f(p)f\left(p\right) 的值为一个关于 pp 的低次多项式,在 pcp^c 处的值 f(pc)f\left(p^c\right) 可以快速计算。

记号

  • P\mathbb P 表示全体质数的集合。
  • 本文一般默认 pPp \in \mathbb P
  • pkp_k 为第 kk 小的质数。特别地,p0=1p_0 = 1
  • lpf(n)\text{lpf}\left(n\right) 表示 nn 的最小质因子(即 Least Prime Factor
  • FP=2pnf(p)F_{\mathbb P} = \sum\limits_{2 \le p \le n} f\left(p\right)
  • Fk=i=2n[pklpf(i)]f(i)F_k = \sum\limits_{i = 2}^n [p_k \le \text{lpf}\left(i\right)] f\left(i\right)

怎么筛?

假设我们已经算出了所有 FF,可以发现答案就是 F1(n)+f(1)=F1(n)+1F_1\left(n\right) + f\left(1\right) = F_1\left(n\right) + 1

求 F

枚举最小质因子及其指数可得

看起来非常的迷,接受不了。
第一个等号后,第一个和式的意思就是枚举最小质因子,并且显然每个合数的最小质因子一定在 n\sqrt n 内。
第二项是补回质数的贡献,因为 c=1c=1 时没有计 f(pc)=f(p)f(p^c)=f(p) 的贡献。

第二个等号后只是把第二项根据定义拆成前缀和的形式。

最后一个等号比较奇怪。
注意到若 picn<pic+1p_i^c \le n < p_i^{c+1},显然有 npic<pi+1\left\lfloor\frac n{p_i^c}\right\rfloor < p_{i+1},故 Fi+1(npic)=0F_{i+1}\left(\left\lfloor\frac n{p_i^c}\right\rfloor\right) = 0

假设已经求得了所有 FPF_{\mathbb P},现在理论上有两种方法计算 FF,根据上式直接递归就是其中一种。

Yet Another F

话说其实 Min_25 筛有无数种实现方式……重要的是埃氏筛的思想。
所以我也不知道从哪里看来的改成递推形式要换 FF 的定义(逃
(如果不换的话也从上式能导出比较显然的递推式,不过也许这种更好写?

Fk(n)=i=2n[pklpf(i)iP]f(i)F_k(n) = \sum\limits_{i=2}^n [p_k \le \text{lpf}(i) \lor i \in \mathbb P]f(i)

易知答案同样为 F1(n)+1F_1(n)+1

考虑一个递推式

于是同样地枚举最小质因子即可。
边界条件为 Fπ(n)+1(n)=FP(n)F_{\pi(\sqrt n)+1}(n) = F_{\mathbb P}(n)
似乎因为指数不会很大所以复杂度仍然是 O(n3/4logn)O\left(\frac{n^{3/4}}{\log n}\right)

求 f 在质数处的值

根据递推式可以发现实际上需要用到的 FPF_{\mathbb P} 只有形如 O(n)O\left(\sqrt n\right) 处。
根据条件,f(p)f\left(p\right) 可以写作 aipi\sum a_i p^i
故我们可以分别讨论每一项的贡献。

Gk,s(n)=i=1n[pk<lpf(i)iP]isG_{k,s}\left(n\right) = \sum\limits_{i = 1}^n [p_k < \text{lpf}\left(i\right) \lor i \in \mathbb P]i^s
熟悉埃筛的话,你会发现这相当于是埃筛筛了 kk 轮后剩下的数的 ss 次方之和。

最后我们想要的实际上就是 Gπ(n),s(n)G_{\pi(\sqrt n),s}\left(n\right),其中 π(n)\pi(\sqrt n) 表示 {p:pn}|\{p : p \le \sqrt n\}|
因为一个合数的最小质因子肯定不会超过 n\sqrt n

考虑转移。
对于 n<pk2n < p_k^2,并没有筛掉新的数,故直接转移即可。
对于 pk2np_k^2 \le n,筛掉的数显然有质因子 pkp_k,故减去 pksGk1,s(npk)p_k^s G_{k - 1,s}\left(\left\lfloor\frac n{p_k}\right\rfloor\right)
然鹅我们筛掉的应该是 lpf(m)=pk\text{lpf}\left(m\right) = p_kmm,故要考虑 lpf(m)<pk\text{lpf}\left(m\right) < p_kmm 的贡献。
实际上就是加回 pksGk1,s(pk1)p_k^s G_{k - 1,s}\left(p_{k - 1}\right)
综上所述,有

复杂度?

并不会证明,由论文可得求 FkF_k 的复杂度为 Θ(n1ϵ)\Theta\left(n^{1-\epsilon}\right)(第一种方法)。
FPF_{\mathbb P} 的复杂度为 O(n3/4logn)O\left(\frac{n^{3/4}}{\log n}\right)(第二种方法的复杂度也长这样)。

一些提示

很多东西可以直接欧拉筛预处理,至于是啥,自己想想。

例题 - LibreOJ 6053.简单的函数

注意到质数除了 22 都是奇数,故 f(p)=p1+2[p=2]f\left(p\right) = p - 1 + 2[p=2]

于是有 G0,0(n)=n+1,G0,1=(n+2)(n1)2G_{0,0}\left(n\right) = -n + 1,G_{0,1} = \dfrac{\left(n+2\right)\left(n-1\right)}2
然后对 22 讨论一下即可。

递归写法的代码见 https://www.alpha1022.me/articles/loj-6053.htm
Yet Another F 中的写法代码:

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
#include <cstdio>
#include <cstring>
#include <cmath>
using namespace std;
const long long N = 1e10;
const int MX = 1e5;
const int mod = 1e9 + 7;
const int inv = 5e8 + 4;
long long n;
int lim;
int vis[MX + 5],prime[MX + 5],cnt,Gprime[MX + 5];
int tot,le[MX + 5],ge[MX + 5];
long long lis[2 * MX + 5];
int G[2 * MX + 5][2],F[2 * MX + 5];
inline int &id(long long x)
{
return x <= lim ? le[x] : ge[n / x];
}
int main()
{
for(register int i = 2;i <= MX;++i)
{
if(!vis[i])
prime[++cnt] = i,Gprime[cnt] = (Gprime[cnt - 1] + i) % mod;
for(register int j = 1;j <= cnt && i * prime[j] <= MX;++j)
{
vis[i * prime[j]] = 1;
if(!(i % prime[j]))
break;
}
}
scanf("%lld",&n),lim = sqrt(n);
for(register long long l = 1,r;l <= n;l = r + 1)
{
r = n / (n / l);
lis[id(n / l) = ++tot] = n / l;
G[tot][0] = (n / l % mod - 1 + mod) % mod,G[tot][1] = (n / l % mod + 2) * (n / l % mod - 1 + mod) % mod * inv % mod;
}
for(register int k = 1;k <= cnt;++k)
{
int p = prime[k];
long long s = (long long)p * p;
for(register int i = 1;lis[i] >= s;++i)
G[i][0] = (G[i][0] - (G[id(lis[i] / p)][0] - (k - 1) + mod) % mod + mod) % mod,
G[i][1] = (G[i][1] - (long long)p * ((G[id(lis[i] / p)][1] - Gprime[k - 1] + mod) % mod) % mod + mod) % mod;
}
for(register int i = 1;i <= tot;++i)
F[i] = (G[i][1] - G[i][0] + mod) % mod;
for(register int k = cnt;k;--k)
{
int p = prime[k];
long long s = (long long)p * p;
for(register int i = 1;lis[i] >= s;++i)
{
long long pw = p,pw2 = s;
for(register int c = 1;pw2 <= lis[i];++c,pw = pw2,pw2 *= p)
F[i] = (F[i] + (long long)(prime[k] ^ c) * (F[id(lis[i] / pw)] - (Gprime[k] - k + mod) + mod) % mod + (prime[k] ^ c + 1)) % mod;
}
}
printf("%d\n",(F[1] + 1 + (n >= 2) * 2) % mod);
}