学习笔记：Min_25筛
Min_25 筛 学习笔记
叫这个名字是因为这是 \(Min\_25\) 神犇发明的，主要用到的思想是对于质数和合数分开计算。
适用范围
对于一个积性函数 \(f(x)\) 求解:
要求：\(f(p)\) ( \(p\) 为质数)是一个关于 \(p\) 的项数较少的多项式或可以快速求值，且 \(f(p^e)\) 可以快速求值
原理
首先考虑求解一个函数 \(G(n)=\sum_{i=1}^{n}[i \in prime]f(i)\)，即 \(n\) 以内所有质数的函数值。
考虑 \(dp\) 求解，定义 \(g(i,m)=\sum_{i=1}^{m}[i \in prime \; or \; minp(i)>p[j]]f'(i)\)
其中 \(minp(i)\) 表示 \(i\) 中最小质因子，\(p[i]\)表示第 \(i\) 个质数。
\(f'(i)\)是一个完全积性函数，满足当 \(i\) 为质数是取值与 \(f(i)\) 相同且前缀和可以快速求解
\(dp\) 边界 \(g(0,m)=\sum_{i=2}^{m}f'(i)\)
当 \(p[i]^2>m\) 时，\(m\) 内不存在最小质因子为 \(p[i]\) 的合数，不会造成贡献，所以 \(g(j,m)=g(j-1,m)\)
当 \(p[i]^2\le m\) 时：
当 \(i-1\) 变成 \(i\) 时，所有最小质因子为 \(p[i]\) 的都被筛掉，因为 \(f'(i)\) 是完全积性函数,提出一个 \(f(p[i])\),剩下的要筛掉的就是最小质因子大于等于\(p[i]\) 的数，因为多算了 \(g(i-1,\lfloor \frac{m}{p[i]} \rfloor )\)里的质数所以要减掉
则 \(G(m)=g(x,m)\)，\(x\) 是小于等于 \(m\) 的质数个数
设 \(S(n,j)=\sum_{i=1}^{n}[minp(i) \ge p[j]]f(i)\)
质数的贡献为 \(G(n)-\sum_{i=1}^{j-1}f(p[i])\)
合数的贡献可以枚举最小质因子和它的次数
实现
实现方面有几个技巧：
1.下标的离散化，因为要用到的 \(\sqrt n\) 个数值域范围很大，所以可以开两个数组 \(id1\) 和 \(id2\)，当值 \(w\) 小于 \(\sqrt n\) 时，下标存在 \(id1[w]\) 里，反之讲下标存在 \(id2[n/w]\) 里，这样两个数组都开成 \(\sqrt n\) 大小就好了
2.求解 \(g(i,m)\) 时可以滚掉第一维
练习
洛谷模板题
将 \(f(x)\) 拆成 \(n^2-n\) 的形式再分别求和，注意到 \(f(x)=x^2\) 和 \(f(x)=x\) 都是完全积性函数，直接Min_25筛求和即可，最后求答案把两个函数减一下
代码
#include <bits/stdc++.h>
#define int long long
#define endl putchar('\n')
using namespace std;
const int M = 1e6+10;
const int inf = 2147483647;
const int Mod = 1e9+7;
const int inv6 = 166666668;
typedef long long ll;
inline int read(){
int x=0,f=0;char c=getchar();
while(!isdigit(c)){
if(c=='-') f=1;
}
do{
x=(x<<1)+(x<<3)+(c^48);
}while(isdigit(c=getchar()));
return f?-x:x;
}
ll p[M],isp[M],id1[M],id2[M],n,tot,sg,cnt;
ll sum1[M],sum2[M],g1[M],g2[M],val[M];
void init(int n){
for(int i=2;i<=n;i++){
if(!isp[i]){
p[++tot]=i;
sum1[tot]=i;sum2[tot]=1ll*i*i%Mod;
sum1[tot]=(sum1[tot]+sum1[tot-1])%Mod;
sum2[tot]=(sum2[tot]+sum2[tot-1])%Mod;
}
for(int j=1;j<=tot&&i*p[j]<=n;j++){
isp[i*p[j]]=1;
if(i%p[j]==0) break;
}
}
}
inline int getid(ll x){
if(x<=sg) return id1[x];
else return id2[n/x];
}
ll s(ll x,ll y){
if(x<=1||x<p[y]) return 0;
int k=getid(x);
ll res=((g2[k]-g1[k]+Mod)%Mod-(sum2[y-1]-sum1[y-1]+Mod)%Mod+Mod)%Mod;
for(int i=y;i<=tot&&p[i]*p[i]<=x;i++){
ll t1=p[i],t2=p[i]*p[i];
for(int j=1;t2<=x;j++,t1=t2,t2*=p[i]){
ll tt1=t1%Mod,tt2=t2%Mod;
res=(res+tt1*(tt1-1)%Mod*s(x/t1,i+1)%Mod+tt2*(tt2-1)%Mod)%Mod;
}
}
return (res+Mod)%Mod;
}
signed main(){
// freopen("a.in","r",stdin);
// freopen("a.out","w",stdout);
n=read();
init(sg=sqrt(n));
for(ll l=1,r;l<=n;l=r+1){
r=n/(n/l);
ll w=n/l;
val[++cnt]=w;
if(w<=sg) id1[w]=cnt;
else id2[n/w]=cnt;
w%=Mod;
g1[cnt]=w*(w+1)/2%Mod;
g2[cnt]=w*(w+1)%Mod*(2ll*w+1ll)%Mod*inv6%Mod;
g1[cnt]=(g1[cnt]-1+Mod)%Mod;
g2[cnt]=(g2[cnt]-1+Mod)%Mod;
}
for(int j=1;j<=tot;j++){
for(int i=1;i<=cnt&&p[j]*p[j]<=val[i];i++){
int k=getid(val[i]/p[j]);
g1[i]=(g1[i]-p[j]*(g1[k]-sum1[j-1]+Mod)%Mod+Mod)%Mod;
g2[i]=(g2[i]-p[j]*p[j]%Mod*(g2[k]-sum2[j-1]+Mod)%Mod+Mod)%Mod;
}
}
printf("%lld\n",(s(n,1)+1)%Mod);
return 0;
}
区间素数个数
求$1 $ ~ $ n$之间素数个数
定义函数 \(f(x)=1\) 和 \(g(x)=1 (x \in prime)\)
发现 \(f(x)\) 是完全积性函数并且在质数处取值与 \(g(x)\) 相同，所以直接用Min_25筛第一部分求 \(g(n)\) 即可
代码
#include <bits/stdc++.h>
#define int long long
#define endl putchar('\n')
using namespace std;
const int M = 1e6+10;
const int inf = 2147483647;
const int Mod = 1e9+7;
const int inv6 = 166666668;
typedef long long ll;
inline int read(){
int x=0,f=0;char c=getchar();
while(!isdigit(c)){
if(c=='-') f=1;
}
do{
x=(x<<1)+(x<<3)+(c^48);
}while(isdigit(c=getchar()));
return f?-x:x;
}
ll p[M],isp[M],val[M],id1[M],id2[M],g[M],n,tot,cnt,sg;
void init(int n){
for(int i=2;i<=n;i++){
if(!isp[i]){
p[++tot]=i;
}
for(int j=1;j<=tot&&i*p[j]<=n;j++){
isp[i*p[j]]=1;
if(i%p[j]==0) break;
}
}
}
inline int getid(int x){
if(x<=sg) return id1[x];
else return id2[n/x];
}
signed main(){
// freopen("a.in","r",stdin);
// freopen("a.out","w",stdout);
n=read();
init(sg=sqrt(n));
for(int l=1,r;l<=n;l=r+1){
r=n/(n/l);
int w=n/l;
val[++cnt]=w;
if(w<=sg) id1[w]=cnt;
else id2[n/w]=cnt;
g[cnt]=w-1;
}
for(int j=1;j<=tot;j++){
for(int i=1;i<=cnt&&p[j]*p[j]<=val[i];i++){
int k=getid(val[i]/p[j]);
g[i]=(g[i]-(g[k]-(j-1)));
}
}
printf("%lld\n",g[getid(n)]);
return 0;
}
某个套路求和题
对于函数 \(f(n)\) 有
求
注意到当 \(x\) 有平方因子时 \(f(x)=0\)，\(x\) 为质数时 \(f(x)=-1\),
其他情况可以看做若干个质因数相乘，如果选奇数个质因数相乘 \(\mu(x)=-1\),而这些数的因数中是奇数个质因数相乘的方案数为偶数（总方案数/2），所以\(f(x)=1\)
所以得到
后一部分Min_25筛，前面设 \(p(i)\) 为 \(i\) 的最大平方因子(\(p^2|n\))，有
然后就没了
代码
#include <bits/stdc++.h>
#define int long long
#define endl putchar('\n')
using namespace std;
const int M = 1e6+10;
const int inf = 2147483647;
const int Mod = 998244353;
const int inv6 = 166666668;
typedef long long ll;
inline int read(){
int x=0,f=0;char c=getchar();
while(!isdigit(c)){
if(c=='-') f=1;
}
do{
x=(x<<1)+(x<<3)+(c^48);
}while(isdigit(c=getchar()));
return f?-x:x;
}
ll p[M],id1[M],id2[M],g[M],mu[M],val[M],tot,cnt,n,sg;
bool isp[M];
void init(int n){
mu[1]=1;
for(int i=2;i<=n;i++){
if(!isp[i]){
p[++tot]=i;
mu[i]=-1;
}
for(int j=1;j<=tot&&p[j]*i<=n;j++){
isp[p[j]*i]=1;
if(i%p[j]==0) break;
mu[i*p[j]]=-mu[i];
}
}
}
inline int idx(int x){
if(x<=sg) return id1[x];
else return id2[n/x];
}
signed main(){
n=read();
init(sg=sqrt(n));
for(int l=1,r;l<=n;l=r+1){
r=n/(n/l);
int w=n/l;
val[++cnt]=w;
if(w<=sg) id1[w]=cnt;
else id2[n/w]=cnt;
g[cnt]=w-1;
}
for(int j=1;j<=tot;j++){
for(int i=1;i<=cnt&&p[j]*p[j]<=val[i];i++){
int k=idx(val[i]/p[j]);
g[i]=g[i]-(g[k]-(j-1));
}
}
ll ans=0;
for(int i=1;i*i<=n;i++) ans+=mu[i]*(n/(i*i));
ans-=2*g[idx(n)];
printf("%lld\n",(ans%Mod+Mod)%Mod);
return 0;
}
质数计数
求满足 \(1 < p \le n\)且 \(p\) 的二进制表示最后两位为\(01\)的质数 \(p\) 有多少个
定义函数 \(f(x)\),当 \(x \mod 4 = 1\) 时 \(f(x)=1\) 当 \(x \mod 4 = 3\) 时 \(f(x)=-1\)，其余情况为 \(0\)
答案为所有 \(x \mod 4 = 1 (x \in prime)\)
注意到 \(f(x)\) 是完全积性函数，Min_25筛处理出素数的前缀和后加上一个所有质数个数减1（减掉2）就是2倍的答案，再除以2即可
代码
#include <bits/stdc++.h>
#define int long long
#define endl putchar('\n')
using namespace std;
const int M = 1e6+10;
const int inf = 2147483647;
const int Mod = 1e9+7;
const int inv6 = 166666668;
typedef long long ll;
inline int read(){
int x=0,f=0;char c=getchar();
while(!isdigit(c)){
if(c=='-') f=1;
}
do{
x=(x<<1)+(x<<3)+(c^48);
}while(isdigit(c=getchar()));
return f?-x:x;
}
ll p[M],isp[M],val[M],id1[M],id2[M],g1[M],g2[M],n,tot,cnt,sg;
void init(int n){
for(int i=2;i<=n;i++){
if(!isp[i]){
p[++tot]=i;
}
for(int j=1;j<=tot&&i*p[j]<=n;j++){
isp[i*p[j]]=1;
if(i%p[j]==0) break;
}
}
}
inline int f(int x){
x%=4;
if(x==1) return 1;
else if(x==3) return -1;
return 0;
}
inline int s(int x){
x%=4;
if(x==1||x==2) return 1;
else return 0;
}
inline int getid(int x){
if(x<=sg) return id1[x];
else return id2[n/x];
}
signed main(){
// freopen("a.in","r",stdin);
// freopen("a.out","w",stdout);
n=read();
init(sg=sqrt(n));
for(int l=1,r;l<=n;l=r+1){
r=n/(n/l);
int w=n/l;
val[++cnt]=w;
if(w<=sg) id1[w]=cnt;
else id2[n/w]=cnt;
g1[cnt]=w-1;
g2[cnt]=s(w)-1;
}
for(int j=1;j<=tot;j++){
for(int i=1;i<=cnt&&p[j]*p[j]<=val[i];i++){
int k=getid(val[i]/p[j]);
g1[i]=(g1[i]-(g1[k]-(j-1)));
g2[i]=(g2[i]-f(p[j])*(g2[k]-g2[getid(p[j-1])]));
}
}
n=getid(n);
printf("%lld\n",(g1[n]+g2[n]-1)>>1);
return 0;
}
