【学习笔记】min_25筛

背景

GDCPC2024
出题人：出这道 min25 筛是给大家增加过题数的 [呲牙][大哭][呲牙][大哭]

min25筛是干啥的

快速求一个积性函数 $F(x)$ 的前缀和
这个 $F(x)$ 需要满足：
$F(p)={\sum }_{i=0}{a}_i{p}^i$，也就说在质数处 $F(p)$ 是关于 $p$ 的多项式
并且 $F(p^k)$ 需要能够快速计算
比杜教筛适用范围更广。
时间复杂度 $O(\frac{n^{\frac{3}{4}}}{logn})$

Step 1 分类

先对 $i$ 按质数和合数分类。
$$\sum _{i=1}^{n}F(i)=F(1)+\sum _{pisaprime}^{n}F(p)+\sum _{xisacomposite}^{n}F(x)$$
由积性函数的性质，可以在合数那里提出一个 $F(p^e)$
$$F(1)+\sum _p^{n}F(p)+\sum _{p\le \sqrt{n},p^e\le n}F(p^e)(\sum _{minp>p,i}^{\lfloor \frac{n}{p}\rfloor }F(i))$$

Step2 预处理

令 $g(n,j)$ 为前 $n$ 个 $F(i)$ 经历 $j$ 次埃氏筛，剩下的 $F(i)$ 值之和。最后我们需要用到的值是
$$g(n,|P|)=\sum _p^{n}F(p)$$
由于 $F(p)$ 是关于 $p$ 的多项式，所以我们只需要令 $F(i)=i^k$ 对应多项式中其中一项，然后再把每一项加起来。即：
$$g(n,j)=\sum _{i=2}^n[iisaprimeorminp(i)>p_j]i^k$$
其中 $minp(i)$ 为 $i$ 的最小质因子，$p_j$为第 $j$ 大的质数。
那么，这玩意怎么求呢？注意，这里的 $n$ 很大，不能直接线性筛。

考虑如何从 $g(n,j-1)$ 转移到 $g(n,j)$，也就是，我们要考虑第 $j$ 次埃氏筛筛掉了什么。
首先，经历了前面 $j-1$ 次埃氏筛之后，剩下的数要么是质数，要么 $minp>p_{j-1}$。所以第 $j$ 次筛掉的数一定是 $minp=p_j$ 的合数。
并且，$p_j$ 能筛掉的最小的数是 $p_j^2$。
我们可以从这些合数提出一个 $p_j$，根据完全积性函数的性质，得到转移式：
g ( n , j ) = { g ( n , j − 1 ) − p j k ( g ( n p j , j − 1 ) − g ( p j − 1 , j − 1 ) ) , p j 2 ≤ n g ( n , j − 1 ) , p j 2 > n g(n,j)=

$$\begin{cases} g(n,j-1)-p_j^k(g(\frac{n}{p_j},j-1)-g(p_{j-1},j-1)), & p_j^2 \leq n \\ g(n,j-1), & p_j^2>n \end{cases}$$ g(n,j)={g(n,j−1)−pjk​(g(pj​n​,j−1)−g(pj−1​,j−1)),g(n,j−1),​pj2​≤npj2​>n​
为什么要减去 $g(p_{j-1},j-1)$?
$g(p_{j-1},j-1)={\sum }_{i=1}^{j-1}{p}_j$，前面的 $g$ 经历 $j-1$ 次埃氏筛会把质数保留下来，而我们只需要减掉合数，所以我们要把质数去掉。
而且，$p_{j-1}$ 是 $\sqrt{n}$ 以内的数，可以使用线性筛预处理前缀和，也就是
$$s{p}_j=\sum _{i=1}^{j}f(p_i)=g(p_j,j)$$
那么我们要的东西$g(n,|P|)$ 的质数集合 $P$ 为不超过 $\sqrt{n}$ 的质数集。
但是 $n$ 还是太大了，我们依旧无法一个一个地把 $g$ 求出来。但是我们又会想到另一个重要结论：
$$\lfloor \frac{\lfloor \frac{n}{a}\rfloor }{b}\rfloor =\lfloor \frac{n}{ab}\rfloor$$
也就说，我们需要的数都能被表示成 $x$ 或者 $\lfloor \frac{n}{x}\rfloor$ 的形式。这些数一共只有 $O(\sqrt{n})$ 个。
那么如何访问这 $O(\sqrt{n})$ 个位置呢？
如果偷懒用 $map$ 的话会多一个 $log$，我们可以用两个数组 $id1,id2$ 存储它们的下标。

if(w[tot]<=sqr) 
	id1[w[tot]]=tot;
else 
	id2[n/w[tot]]=tot;
1
2
3
4

step3 求答案

我们令 $S(n,j)$ 为 $2$~$n$ 中 $minp>p_j$ 的函数值之和，即：
$$S(n,j)=\sum _{i=2}^n[minp(i)>p_j]$$
那么我们要求的答案就是 $S(n,0)$

同样的，我们依旧把和分成质数的和，合数的和。枚举它们最小的因子，那么就有下面的转移：
$$S(n,j)=g(n,|P|)-s{p}_j+\sum _{p_k^e\le n\&k>x}f(p_k^e)(S(\frac{n}{p_k^e},k)+[e\ne 1])$$
前面的 $g(n,|P|)-s{p}_j$ 显然是质数的答案部分。
后面的合数同样可以根据积性函数的性质提出 $f(p_k^e)$
需要注意的是，$S(\frac{n}{p_k^e},k)$ 是会把 $p_k^e[e>1]$ 筛掉的，所以要把它给补上，也就是加上 $[e\ne 1]$

至此，min25筛的流程就讲完啦。

例题

【模板】Min_25 筛

思路

显然，$f(p)=p^2-p$，那我们分别对它的二次项，一次项进行预处理，然后求解即可。
记得最后把 $f(1)$ 加上

代码

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=2e5+7,inf=1e18,mod=1e9+7;
vector<int> g1(N),g2(N),id1(N),id2(N),sp1(N),sp2(N),p,w(N);
int sqr,tot=0,n;
int power(int x,int t)
{
	int b=1;
	while(t)
	{
		if(t&1) b=b*x%mod;
		x=x*x%mod; t>>=1;
	}
	return b;
}
void init(int n)
{
	vector<bool> bz(n+1);
	p.push_back(0);
	tot=0;
	bz[1]=1;
	for(int i=1; i<=n; i++)
	{
		if(!bz[i])
		{
			p.push_back(i);
			int now=p.size()-1;
			sp1[now]=(sp1[now-1]+i)%mod;
			sp2[now]=(sp2[now-1]+i*i%mod)%mod;
		}
		for(auto j:p)
		{
			if(!j) continue;
			if(i*j>n) break;
			bz[i*j]=1;
			if(i%j==0) break;
		}
	}
}
int S(int x,int y)
{
	if(x<=1||p[y]>=x) return 0;
	int k=(x<=sqr)?id1[x]:id2[n/x];
	//prime part
	int ans=(g2[k]-g1[k])-(sp2[y]-sp1[y]);
	ans%=mod;
	//other parts
	for(int k=y+1; k<p.size()&&p[k]*p[k]<=x; k++)
	{
		int t=p[k];
		for(int e=1; t<=x; e++,t*=p[k])
		{
			int p=t%mod;// avoid overflow
			(ans+=p*(p-1)%mod*(S(x/t,k)+(e!=1))%mod)%=mod;
		}
	}
	return ans;
}
void O_o()
{
	//f(p^k) = p^k(p^k - 1)
	cin>>n;
	sqr=sqrt(n);
	init(sqr);
	//let g(n,0) = \sum_{i=1}^n f'(i) 
	int inv2=power(2,mod-2),inv3=power(3,mod-2);
	for(int i=1,j; i<=n; i=j+1)
	{
		j=n/(n/i);
		w[++tot]=n/i;
		int now=w[tot]%mod;
		g1[tot]=now*(now+1)%mod*inv2%mod-1;
		g2[tot]=now*(now+1)%mod*(2*now+1)%mod*inv2%mod*inv3%mod-1;// 4+9+16...
		if(w[tot]<=sqr) id1[w[tot]]=tot;
		else id2[n/w[tot]]=tot;
	}
	//get g(n,|P|)
	for(int i=1; i<p.size(); i++)
	{
		for(int j=1; j<=tot&&p[i]*p[i]<=w[j]; j++)
		{
			//get g(w[j],i)
			int k=w[j]/p[i]<=sqr?id1[w[j]/p[i]]:id2[n/(w[j]/p[i])];
			//g(w[j],i) = g(w[j],i-1) - f'(p[i])*(g(w[k],j-1)-sp[i-1])
			(g1[j]-=p[i]*(g1[k]-sp1[i-1]+mod)%mod)%=mod;
			(g2[j]-=p[i]*p[i]%mod*(g2[k]-sp2[i-1]+mod)%mod)%=mod;
		}
	}
	int ans=S(n,0)+1;
	(ans+=mod)%=mod;
	cout<<ans<<"\n";
}
signed main()
{
	ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
	cout<<fixed<<setprecision(2);
	int T=1;
//	cin>>T;
	while(T--)
	{
		O_o();
	}
}
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
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104

简单的函数

思路

先考虑质数处的值，显然除了 2 以外，所有质数都是奇数。但我们可以先令 $f(p)=p-1$，最后答案再加上 2。
那么我们需要对 $p$ 和 $-1$ 进行预处理，然后直接套 min25 筛即可。

代码

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=2e5+7,inf=1e18,mod=1e9+7;
vector<int> g1(N),g0(N),sp1(N),id1(N),id2(N),p,w(N);
int n,sqr,tot;
void init(int n)
{
	p.push_back(0);
	tot=0;
	vector<bool> bz(n+1);
	for(int i=2; i<=n; i++)
	{
		if(!bz[i])
		{
			p.push_back(i);
			int now=p.size()-1;
			sp1[now]=(sp1[now-1]+i)%mod;
		}
		for(auto j:p)
		{
			if(!j) continue;
			if(i*j>n) break;
			bz[i*j]=1;
			if(i%j==0) break;
		}
	}
}
int S(int x,int y)
{
	if(x<=1||p[y]>=x) return 0;
	int k=(x<=sqr)?id1[x]:id2[n/x];
	int ans=g1[k]-g0[k]-sp1[y]+y;
	ans%=mod;
	for(int k=y+1; k<p.size()&&p[k]*p[k]<=x; k++)
	{
		int t=p[k];
		for(int e=1; t<=x; e++,t*=p[k])
		{
			(ans+=(p[k]^e)*(S(x/t,k)+(e!=1))%mod)%=mod;
		}
	}
	return ans;
}
void O_o()
{
	cin>>n;
	if(n==1)
	{
		cout<<1<<"\n";
		return;
	}
	sqr=sqrt(n);
	init(sqr);
	for(int i=1,j; i<=n; i=j+1)
	{
		j=n/(n/i);
		w[++tot]=n/i;
		int now=w[tot]%mod;
		g1[tot]=now*(now+1)/2%mod-1;
		g0[tot]=now-1;
		if(w[tot]<=sqr) id1[w[tot]]=tot;
		else id2[n/w[tot]]=tot;
	}
	for(int i=1; i<p.size(); i++)
	{
		for(int j=1; j<=tot,p[i]*p[i]<=w[j]; j++)
		{
			int k=w[j]/p[i]<=sqr?id1[w[j]/p[i]]:id2[n/(w[j]/p[i])];
			//g(w[j],i) = g(w[j],i-1) - f'(p[i])*(g(w[k],j-1)-sp[i-1])
			(g1[j]-=p[i]*(g1[k]-sp1[i-1])%mod)%=mod;
			(g0[j]-=(g0[k]-(i-1))%mod)%=mod;
		}
	}
	int ans=S(n,0)+3;
	(ans+=mod)%=mod;
	cout<<ans<<"\n";
}
signed main()
{
	ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
	cout<<fixed<<setprecision(2);
	int T=1;
//	cin>>T;
	while(T--)
	{
		O_o();
	}
}
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
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89

【GDCPC2024】J.另一个计数问题

题解戳这里喵~