专栏文章

Min_25 筛算法学习笔记

HHPXXZYY
算法·理论参与者 1已保存评论 0

文章操作

快速查看文章及其快照的属性，并进行相关操作。

当前评论
0 条
当前快照
1 份
快照标识符
@miots7x6
此快照首次捕获于
2025/12/03 01:01
3 个月前
此快照最后确认于
2025/12/03 01:01
3 个月前
查看原文
牛客暑期多校连续两场都考了 Min_25 筛，看来不得不学了。

前言

Min_25 筛是一种计算积性函数前缀和的高效算法。
时间复杂度好像是 O(n34log⁡n)O\left ( \frac{n^{\frac{3}{4}}}{\log n} \right ) 的。

模板题（洛谷 P5325）

[Problem Discription]\color{blue}{\texttt{[Problem Discription]}}
已知某积性函数 f(x)f(x) 满足 f(pk)=pk(pk−1)f(p^{k})=p^{k}(p^{k}-1)，其中 pp 为质数，k∈N+k \in \mathbb{N}_{+}。求 ∑i=1nf(i)\sum\limits_{i=1}^{n} f(i)。
1≤n≤1×10101 \leq n \leq 1 \times 10^{10}。
[Analysis]\color{blue}{\texttt{[Analysis]}}
首先，我们先化简问题，考虑求幂函数 F(p)=pkF(p)=p^{k} 的前缀和。
注，对于幂函数而言，F(x)F(x) 不仅是积性函数，更是完全积性函数。想一想，怎么证明。
先把求和式子按质数、合数分类：
∑i=1nF(i)=∑p is prime, p≤nF(p)+∑i=1n[i is not a prime]F(i)\sum\limits_{i=1}^{n} F(i)=\sum\limits_{p \text{ is prime, } p \leq n}F(p)+\sum\limits_{i=1}^{n}[i \text{ is not a prime}]F(i)
然后枚举最小质因子：
∑i=1nF(i)=∑p is prime, p≤nF(p)+∑1≤p≤n,1≤pe≤nF(pe)(∑i=1⌊npe⌋[minp(i)≥p]F(i))\sum\limits_{i=1}^{n} F(i)=\sum\limits_{p \text{ is prime, } p \leq n}F(p)+\sum\limits_{1 \leq p \leq \sqrt{n},1 \leq p^{e}\leq n}F(p^{e})\left ( \sum\limits_{i=1}^{\left \lfloor \frac{n}{p^{e}} \right \rfloor} [\text{minp}(i) \geq p]F(i) \right )
minp(i)\text{minp}(i) 表示 ii 的最小质因子，这里利用了积性函数的性质。
考虑一个天才的 dp。记 pip_{i} 表示从小到大的第 ii 个素数， g(n,j)g(n,j) 表示前 nn 个自然数中满足本身是素数或者最小质因数大于 pjp_{j} 的所有整数 ii 的 kk 次幂的和，即：
g(n,j)=∑1≤i≤nik[i∈Prime or minp(i)>pj]g(n,j)=\sum\limits_{1 \leq i \leq n}i^{k}\left [ i \in \mathbb{P}\text{rime or minp}(i) > p_{j} \right ]
从 g(n,j−1)g(n,j-1) 转移到 g(n,j)g(n,j) 时，有一些整数不满足条件，需要被剔除。容易发现，被剔除的整数满足最小质因数恰好为 pjp_{j}。
可以提出一个公因子 pjp_{j}，剩下的部分就不能有比 pjp_{j} 更小的质因子了，即 g(⌊npj⌋,j)−g(pj−1,j−1)g\left ( \left \lfloor \frac{n}{p_{j}} \right \rfloor,j \right )-g(p_{j-1},j-1)，后面那一项是为了把质数全部删除。
这样就得到了 gg 的递推关系式：
g(n,j)=g(n,j−1)−pjk[g(⌊npj⌋,j)−g(pj−1,j−1)]g(n,j)=g(n,j-1)-p_{j}^{k} \left [ g\left ( \left \lfloor \frac{n}{p_{j}} \right \rfloor,j \right )-g(p_{j-1},j-1) \right ]
这里体现出来了完全积性函数的好处，虽然只提取了一个 pjp_{j}，但是函数值可以直接相乘而不需要考虑是否互质。
注意到 g(pj−1,j−1)g(p_{j-1},j-1) 其实就是前 (j−1)(j-1) 个质数的 kk 次幂的和，且 pj−1≤np_{j-1} \leq \sqrt{n}，直接线性筛即可，对应于代码中的 pre 数组。
11 到 nn 中所有质数的 kk 次幂的和其实就是 g(n,x)g(n,x)，其中 pxp_{x} 是最后一个小于等于 n\sqrt{n} 的质数。为了方便，把 g(n,x)g(n,x) 简单记为 g(n)g(n)。
为了求出 g(n)g(n)，我们需要用到以下性质：
⌊⌊na⌋b⌋=⌊nab⌋\left \lfloor \frac{\left \lfloor \frac{n}{a} \right \rfloor}{b} \right \rfloor=\left \lfloor \frac{n}{ab} \right \rfloor
因此，无论做了多少次除法，最终得到的数一定可以写成 ⌊nx⌋\left \lfloor \frac{n}{x} \right \rfloor 的形式。可以除法分块。
为了把这 O(n)O(\sqrt{n}) 个数保存下来，我们需要开两个下标数组，具体看代码。
接下来就是求解原问题了。
还是考虑相似的 dp。答案可以分成质数的函数和和合数的函数和两部分。质数的函数和可以用多项式直接相加减得到。考虑合数的函数和。
用 S(n,j)S(n,j) 表示
S(n,j)=∑i=1nf(i)[minp(i)>pj]S(n,j)=\sum\limits_{i=1}^{n}f(i)[\text{minp}(i)>p_{j}]
则
S(n,j)=g(n)−pre(x)+∑pke≤n,k>xf(pkx){S(npke,k)+[e≠1]}S(n,j)=\color{blue}{g(n)-\text{pre}(x)}\text{\color{black}{+}}\color{red}{\sum\limits_{p_{k}^{e}\leq n,k>x}f(p_{k}^{x}) \left \{ S \left ( \frac{n}{p_{k}^{e}},k \right ) + [e\not = 1] \right \}}
蓝色项是大于 pjp_{j} 的质数的贡献，红色项是最小质因数大于 pjp_{j} 的合数的贡献。
具体参看模板题的代码实现。
Code\color{blue}{\text{Code}}
CPP
int prime[N],prcnt;
bool is_prime[N];
ll pre1[N],pre2[N];

void get_prime(int n){
	fill(is_prime+2,is_prime+n+1,true);
	
	for(int i=2;i<=n;i++){
		if (is_prime[i]){
			prime[++prcnt]=i;
			pre1[prcnt]=(pre1[prcnt-1]+i)%mod;
			pre2[prcnt]=(pre2[prcnt-1]+1ll*i*i%mod)%mod;
		}
		
		for(int j=1;j<=prcnt;j++){
			if (1ll*i*prime[j]>n) break;
			is_prime[i*prime[j]]=false;
			if (i%prime[j]==0) break;
		}
	}
}

ll g1[N<<1],g2[N<<1],w[N<<1];
int ind1[N],ind2[N];

ll n;int m,wcnt;

const int inv2=500000004;
const int inv3=333333336;

void get_g_and_ind(){
	for(ll l=1;l<=n;){
		ll r=n/(n/l);
		
		w[++wcnt]=n/l;
		
		int tmp=(n/l)%mod;
		g1[wcnt]=1ll*tmp*(tmp+1)%mod*inv2%mod-1;
		g2[wcnt]=1ll*tmp*(tmp+1)%mod*(2*tmp+1)%mod*inv2%mod*inv3%mod-1;
		
		if (n/l<=m) ind1[n/l]=wcnt;
		else ind2[r]=wcnt;
		
		l=r+1;
	}
	
	for(int i=1;i<=prcnt;i++)
		for(int j=1;j<=wcnt;j++){
			if (1ll*prime[i]*prime[i]>w[j]) break;
			
			int k=(w[j]/prime[i]<=m)?ind1[w[j]/prime[i]]:ind2[n/(w[j]/prime[i])];
			
			g1[j]=(g1[j]-1ll*prime[i]*(g1[k]-pre1[i-1]+mod)%mod+mod)%mod;
			g2[j]=(g2[j]-1ll*prime[i]*prime[i]%mod*(g2[k]-pre2[i-1]+mod)%mod+mod)%mod;
		}//这里 g 是可以滚动数组的，最终保留的信息就是 g(n,x)
}

int S(ll x,int y){
	if (prime[y]>=x) return 0;
	
	int k=(x<=m)?ind1[x]:ind2[n/x];
	int ans=((g2[k]-g1[k]+mod)%mod-(pre2[y]-pre1[y]+mod)%mod+mod)%mod;
	
	for(int i=y+1;i<=prcnt;i++){
		if (1ll*prime[i]*prime[i]>x) break;
		
		ll pe=prime[i];
		for(int e=1;pe<=x;e++,pe*=prime[i]){
			int tmp=pe%mod;
			ans=(ans+1ll*tmp*(tmp-1)%mod*(S(x/pe,i)+(e!=1))%mod)%mod;
		}
	}
	
	return ans;
}

int main(){
	cin>>n;
	m=sqrt(n);
	
	get_prime(m);
	get_g_and_ind();
	
	cout<<(S(n,0)+1)%mod;
	
	return 0;
}

本质理解

虽然 Min_25 筛正宗的用法是用来求积性函数的，但是有些非积性函数的问题也可以用 Min_25 筛来做。
由上面的代码我们可以知道，我们其实不关心 g(n,j)g(n,j) 的值是多少，我们关心的是 g(n,x)g(n,x)，即 g(n)g(n)。而 g(n)g(n) 里保存的是所有的质数的信息，而且与合数完全无关。
而在 SS 的转移式中，蓝色部分求出来的也是函数在所有质数点的取值的和，只有红色部分包含合数，且红色部分只包含合数。
于是，在某些特定的问题中，我们只关心函数在质数的和，或函数在合数的和（两者等价），我们一样可以用 Min_25 筛解决。
一个经典的问题是求某个区间 [l,r][l,r] 中质数的个数，可以用前缀和的方法转化为求 [1,n][1,n] 中质数的个数。
定义函数
1 \quad n \text{ is a Prime}\\ 0 \quad n \text{ isn't a Prime} \end{cases}$$ 显然这个函数不是积性函数。 但是我们发现 $F(n)=1$ 是积性函数，那么我们完全可以就求 $F(n)$ 的前缀和，不过是求在质数点的函数的和。也就是只要蓝色部分即可。 ____________________ $2025$ 牛客暑期多校 $\text{Con } 2 \text{ E}$ 题。 $\color{blue}{\texttt{[Problem Discription]}}$ 求 $[l,r]$ 中有多少个数可以加上自己的一个质因数得到一个完全平方数。 $2 \leq l \leq r \leq 1 \times 10^{18}$。 $\color{blue}{\texttt{[Analysis]}}$ 可以证明，满足条件的数一定可以写成如下形式： $$x^{2}-p$$ 其中 $p$ 是 $x$ 的一个质因数，且可以证明不同的 $x$ 生成的一定是不同的整数。 转化为求 $[1,n]$。$\left \lfloor \sqrt{n} \right \rfloor+1$ 是可以生成合法的解的，但是就一个数，特判就好。 对于 $1 \dots \sqrt{n}$ 中的数，其实就是求质因数的个数的和，即 $$\sum\limits_{1 \leq p \leq \sqrt{n}} \left \lfloor \frac{n}{p} \right \rfloor$$ 对于合适范围内的 $p$，直接这么干就好了。对于很大的 $p$，$\left \lfloor \frac{n}{p} \right \rfloor$ 很小，考虑枚举 $\left \lfloor \frac{n}{p} \right \rfloor$，得到一个区间，求这个区间的质数个数就好。 $\color{blue}{\text{Code}}$ ```cpp const int N=1e5+100; #define ll long long int prime[N],prcnt,pre[N]; bool is_prime[N]; void get_prime(int n){ fill(is_prime+2,is_prime+n+1,true); for(int i=2;i<=n;i++){ if (is_prime[i]){ prime[++prcnt]=i; pre[prcnt]=prcnt; } for(int j=1;j<=prcnt;j++){ if (1ll*prime[j]*i>n) break; is_prime[i*prime[j]]=false; if (i%prime[j]==0) break; } } } int w[N<<1],ind1[N],ind2[N],wcnt; ll g[N<<1];int n,m; void get_g_and_ind(){ for(int l=1;l<=n;){ int tmp=n/l,r=n/tmp; w[++wcnt]=tmp; g[wcnt]=tmp-1; if (tmp<=m) ind1[tmp]=wcnt; else ind2[r]=wcnt; l=r+1; } for(int i=1;i<=prcnt;i++) for(int j=1;j<=wcnt;j++){ if (1ll*prime[i]*prime[i]>w[j]) break; int k=(w[j]/prime[i]<=m)?ind1[w[j]/prime[i]]:ind2[n/(w[j]/prime[i])]; g[j]-=g[k]-pre[i-1]; } } ll S(ll x,int y){ if (prime[y]>=x) return 0; int k=(x<=m)?ind1[x]:ind2[n/x]; return g[k]-pre[y]; } void init(){ memset(pre,0,sizeof(pre)); memset(g,0,sizeof(g)); memset(ind1,0,sizeof(ind1)); memset(ind2,0,sizeof(ind2)); wcnt=prcnt=0; } ll Count(int x){ n=x; m=sqrt(n); init(); get_prime(m); get_g_and_ind(); ll res=0; for(int i=1;i<=prcnt;i++) res+=x/prime[i]; for(int i=1;i<=m;i++) if (x/i>prime[prcnt]) res+=S(x/i,0)-prcnt; return res; } ll calc(ll x){ int r=sqrt(x); ll ans=Count(r); int tmp=r+1; for(int i=2;1ll*i*i<=tmp;i++) if (tmp%i==0){ if (1ll*(r+1)*(r+1)-i<=x) ans++; while (tmp%i==0) tmp/=i; } if (tmp>1&&1ll*(r+1)*(r+1)-tmp<=x) ans++; return ans; } int main(){ ll l,r; cin>>l>>r; cout<<calc(r)-calc(l-1); return 0; } ```

评论

共 0 条评论，欢迎与作者交流。

正在加载评论...