之前考试就学了一下 m i n _ 25 min\_25 min_25筛，现在理解得更多一些，补一波笔记
首先要知道这是一种可以筛出积性函数前缀和的东西，大概可以做洲阁筛能做的东西，复杂度是 O ( n 3 4 l o g n ) O(\frac{n^{\frac{3}{4}}}{logn}) O(lognn43)的，大概能做到 1 0 10 10^{10} 1010~ 1 0 11 10^{11} 1011
算法步骤
设
f
(
i
)
f(i)
f(i)是一个积性函数，或者说
f
(
i
)
f(i)
f(i)是一个关于质数
p
p
p的多项式，并且
f
(
p
q
)
f(p^q)
f(pq)可以快速求出来
要求
∑
i
=
1
n
f
(
i
)
\sum_{i=1}^nf(i)
∑i=1nf(i)
（下面内容中 m i n p x minp_x minpx表示 x x x的最小质因子， p j p_j pj表示筛出来的第 j j j个质数
- 线性筛出 n \sqrt{n} n以内的所有质数和 f ( i ) f(i) f(i)
- 求出
∑
i
=
1
n
[
i
∈
p
r
i
m
e
]
f
(
i
)
\sum_{i=1}^n[i\in prime]f(i)
∑i=1n[i∈prime]f(i)
设 g ( n , j ) g(n,j) g(n,j)表示 ∑ i = 1 n [ i ∈ p i r m e o r m i n p i > p j ] f ( i ) \sum_{i=1}^n[i\in pirme\ or\ minp_i>p_j]f(i) ∑i=1n[i∈pirme or minpi>pj]f(i)
考虑从 g ( n , j − 1 ) − > g ( n , j ) g(n,j-1)->g(n,j) g(n,j−1)−>g(n,j)，需要把 p j p_j pj的贡献算出来减掉，就是
g ( n , j ) = { g ( n , j − 1 ) p j 2 > n 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)=\begin{cases}g(n,j-1)&p_j^2>n\\g(n,j-1)-f(p_j)\times [g(\frac{n}{p_j},j-1)-\sum_{i=1}^{j-1}f(p_i)]&p_j^2\le n\end{cases} g(n,j)={g(n,j−1)g(n,j−1)−f(pj)×[g(pjn,j−1)−∑i=1j−1f(pi)]pj2>npj2≤n
上式就是考虑 m i n p = p j minp=p_j minp=pj的那些数的贡献，因为是积性函数，所以可以把 p j p_j pj提出来，剩下有
g ( n p j , j − 1 ) g(\frac{n}{p_j},j-1) g(pjn,j−1)那么多，其中不包括质数的所以要减掉 ∑ i = 1 j − 1 f ( p i ) \sum_{i=1}^{j-1}f(p_i) ∑i=1j−1f(pi) - 求出
∑
i
=
1
n
f
(
i
)
\sum_{i=1}^nf(i)
∑i=1nf(i)
这一步中其实只需要考虑合数的情况，因为质数已经都算出来了
设 S ( n , j ) = ∑ i = 1 n [ m i n p i ≥ p j ] f ( i ) S(n,j)=\sum_{i=1}^n[minp_i\ge p_j]f(i) S(n,j)=∑i=1n[minpi≥pj]f(i)，考虑每个合数在 m i n p > p j minp>p_j minp>pj时的贡献，可以枚举 m i n p minp minp和指数 e e e，就有
S ( n , j ) = g ( n , j ) − ∑ i = 1 j − 1 f ( p i ) + ∑ k = j p k 2 ≤ n ∑ e = 1 p k e + 1 ≤ n S ( n p k e , k + 1 ) × f ( p k e ) + f ( p k e + 1 ) S(n,j)=g(n,j)-\sum_{i=1}^{j-1}f(p_i)+\sum_{k=j}^{p_k^2\le n}\sum_{e=1}^{p_k^{e+1}\le n}S(\frac{n}{p_k^e},k+1)\times f(p_k^e)+f(p_k^{e+1}) S(n,j)=g(n,j)−i=1∑j−1f(pi)+k=j∑pk2≤ne=1∑pke+1≤nS(pken,k+1)×f(pke)+f(pke+1)
答案就是 S ( n , 1 ) + f ( 1 ) = ∑ i = 1 n f ( i ) S(n,1)+f(1)=\sum_{i=1}^nf(i) S(n,1)+f(1)=∑i=1nf(i)
这个通常是递归求解，也有非递归版的，有时会比递归版快很多，下面有
例题
LOJ#6053. 简单的函数
大概是
m
i
n
_
25
min\_25
min_25筛模板题？
发现除了
p
=
2
p=2
p=2是
f
(
p
)
=
3
f(p)=3
f(p)=3，其余情况
f
(
p
)
=
p
−
1
f(p)=p-1
f(p)=p−1，所以算出
g
(
n
)
=
∑
i
=
1
n
[
i
∈
p
r
i
m
e
]
i
g(n)=\sum_{i=1}^n[i\in prime]i
g(n)=∑i=1n[i∈prime]i和
h
(
n
)
=
∑
i
=
1
n
[
i
∈
p
r
i
m
e
]
1
h(n)=\sum_{i=1}^n[i\in prime] 1
h(n)=∑i=1n[i∈prime]1就可以算出所有
f
(
p
)
f(p)
f(p)的前缀和，进而算出
S
(
n
,
1
)
S(n,1)
S(n,1)，另外要特判一下含
2
2
2的情况，要
+
2
+2
+2
代码如下：
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#define N 1000005
#define LL long long
using namespace std;
const int mod=1e9+7;
LL n,w[N];
int cnt,pri[N],sqr,tot,m,id1[N],id2[N],g[N],h[N],sum[N];
bool vis[N];
inline int qpow(int x,int k){
int ret=1;
while(k){
if(k&1) ret=1LL*ret*x%mod;
x=1LL*x*x%mod; k>>=1;
} return ret;
}
inline void get_prime(int n){
for(int i=2;i<=n;i++){
if(!vis[i]) pri[++cnt]=i,sum[cnt]=(sum[cnt-1]+i)%mod;
for(int j=1;j<=cnt && i*pri[j]<=n;j++){
vis[i*pri[j]]=1;
if(i%pri[j]==0) break;
}
}
}
LL S(LL x,int y){
if(x<=1 || pri[y]>x) return 0;
int pos=(x<=sqr)?id1[x]:id2[n/x];
LL ret=g[pos]-h[pos]-sum[y-1]+(y-1); ret=(ret+mod)%mod;
if(y==1) ret+=2;
for(int i=y;i<=cnt && 1LL*pri[i]*pri[i]<=x;i++){
LL p1=pri[i],p2=1LL*pri[i]*pri[i];
for(int e=1;p2<=x;e++,p1=p2,p2*=pri[i])
(ret+=1LL*S(x/p1,i+1)*(pri[i]^e)%mod+(pri[i]^(e+1))%mod)%=mod;
} return ret;
}
int main(){
scanf("%lld",&n); int tmp=qpow(2,mod-2);
sqr=sqrt(n); get_prime(sqr);
for(LL i=1,j;i<=n;i=j+1){
j=n/(n/i); w[++m]=n/i;
if(w[m]<=sqr) id1[w[m]]=m;
else id2[n/w[m]]=m;
h[m]=(w[m]-1)%mod;
g[m]=((w[m]+2)%mod)*((w[m]-1)%mod)%mod*tmp%mod;
}
for(int j=1;j<=cnt;j++)
for(int i=1;i<=m&&1LL*pri[j]*pri[j]<=w[i];i++){
int k=(w[i]/pri[j]<=sqr)?id1[w[i]/pri[j]]:id2[n/(w[i]/pri[j])];
g[i]=(g[i]-1LL*pri[j]*(g[k]-sum[j-1])%mod+mod)%mod;
h[i]=(h[i]-(h[k]-(j-1))+mod)%mod;
}
printf("%lld\n",S(n,1)+1);
return 0;
}
UOJ#188. 【UR #13】Sanrd
f
(
x
)
f(x)
f(x)表示
x
x
x的次大质因子，看起来不是积性函数，并且也不是关于
p
p
p的多项式
但它也是可以用 m i n _ 25 min\_25 min_25筛的！！！
首先可以考虑上面说的
S
S
S函数怎么求，是要枚举
m
i
n
p
=
p
j
,
p
j
e
minp=p_j,p_j^e
minp=pj,pje考虑这些数的贡献，那实际上这个问题也可以这么解决，只考虑质因数
≤
2
\le 2
≤2的时候，因为
>
2
>2
>2可以继续递归求解，如果次小质因子
=
p
j
=p_j
=pj，那么贡献其实就是剩下的数中是
p
r
i
m
e
prime
prime且
≥
p
j
\ge p_j
≥pj的个数，然后就是前面那种套路的想法了
S
(
n
,
j
)
=
∑
k
=
j
p
k
2
≤
n
∑
e
=
1
p
k
e
+
1
≤
n
S
(
n
p
k
e
,
k
+
1
)
+
p
k
×
∑
i
=
p
k
n
p
k
e
[
i
∈
p
r
i
m
e
]
S(n,j)=\sum_{k=j}^{p_k^2\le n}\sum_{e=1}^{p_k^{e+1}\le n}S(\frac{n}{p_k^e},k+1)+p_k\times \sum_{i=p_k}^{\frac{n}{p_k^e}}[i\in prime]
S(n,j)=k=j∑pk2≤ne=1∑pke+1≤nS(pken,k+1)+pk×i=pk∑pken[i∈prime]
代码如下：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#define N 1000005
#define LL long long
using namespace std;
LL l,r,w[N],g[N],n;
int cnt,pri[N],m,id1[N],id2[N],sqr;
bool vis[N];
inline void get_prime(int n){
for(int i=2;i<=n;i++){
if(!vis[i]) pri[++cnt]=i;
for(int j=1;j<=cnt && i*pri[j]<=n;j++){
vis[i*pri[j]]=1;
if(i%pri[j]==0) break;
}
}
}
inline LL S(LL x,int y){
if(x<=2||pri[y]>x) return 0;
LL ret=0;
for(int i=y;i<=cnt && 1LL*pri[i]*pri[i]<=x;i++)
for(LL t=pri[i];t*pri[i]<=x;t*=pri[i]){
int pos=(x/t<=sqr)?id1[x/t]:id2[n/(x/t)];
ret+=S(x/t,i+1)+1LL*pri[i]*(g[pos]-(i-1));
}
return ret;
}
inline LL solve(LL nn){
sqr=sqrt(nn); m=0; n=nn;
for(LL i=1,j;i<=n;i=j+1){
j=n/(n/i); w[++m]=n/i;
if(w[m]<=sqr) id1[w[m]]=m;
else id2[n/w[m]]=m;
g[m]=w[m]-1;
}
for(int i=1;i<=cnt && 1LL*pri[i]*pri[i]<=n;i++)
for(int j=1;j<=m && 1LL*pri[i]*pri[i]<=w[j];j++){
int pos=(w[j]/pri[i]<=sqr)?id1[w[j]/pri[i]]:id2[n/(w[j]/pri[i])];
g[j]-=g[pos]-(i-1);
}
return S(n,1);
}
int main(){
scanf("%lld%lld",&l,&r);
get_prime((int)sqrt(r));
printf("%lld\n",solve(r)-solve(l-1));
return 0;
}
LOJ#572. 「LibreOJ Round #11」Misaka Network 与求和
这道题其实就是上一道题套了个莫比乌斯反演
套路反演得到一个可以数论分块的式子，然后中间那部分再两边卷一下
1
1
1函数，
μ
∗
1
=
e
\mu*1=e
μ∗1=e，右边就没了，然后左边移项就可以得到
g
(
n
)
=
∑
i
=
1
n
F
(
i
)
−
∑
i
=
2
n
g
(
n
i
)
g(n)=\sum_{i=1}^nF(i)-\sum_{i=2}^ng(\frac{n}{i})
g(n)=i=1∑nF(i)−i=2∑ng(in)
发现
F
F
F前缀和就用上一个题的方法
m
i
n
_
25
min\_25
min_25筛筛出来，然后再数论分块一下就可以算这个式子了
递归版复杂度玄学，但可以跑过
非递归版复杂度十分优秀，大概就是
O
(
m
i
n
_
25
筛
)
O(min\_25筛)
O(min_25筛)级别的
代码如下：
（中间注释掉那部分是递归版的，非递归好像能跑递归的
1
3
\frac{1}{3}
31
#include<bits/stdc++.h>
#define N 1000005
#define uint unsigned int
using namespace std;
int K,m,cnt,pri[N],id1[N],id2[N],sqr;
bool vis[N];
uint n,w[N],g[N],f[N],h[N],ans,sum[N];
template<class T>inline void rd(T &x){
x=0; short f=1; char c=getchar();
while(c<'0' || c>'9') f=c=='-'?-1:1,c=getchar();
while(c<='9' && c>='0') x=x*10+c-'0',c=getchar();
x*=f;
}
inline int get(uint x){return x<=sqr?id1[x]:id2[n/x];}
inline uint qpow(uint x,int k){
uint ret=1;
while(k){
if(k&1) ret*=x;
x*=x; k>>=1;
} return ret;
}
inline void get_prime(int n){
for(int i=2;i<=n;i++){
if(!vis[i]) pri[++cnt]=i,f[cnt]=qpow(i,K);
for(int j=1;j<=cnt&&i*pri[j]<=n;j++){
vis[i*pri[j]]=1;
if(i%pri[j]==0) break;
}
}
}
inline void init(){
for(uint i=1,j;i<=n;i=j+1){
j=n/(n/i); w[++m]=n/i;
if(w[m]<=sqr) id1[w[m]]=m;
else id2[n/w[m]]=m;
g[m]=w[m]-1;
}
for(int i=1;i<=cnt && pri[i]<=sqr;i++)
for(int j=1;j<=m&&1LL*pri[i]*pri[i]<=w[j];j++)
g[j]-=g[get(w[j]/pri[i])]-(i-1);
for(int i=cnt;i;i--)
for(int j=1;j<=m && 1LL*pri[i]*pri[i]<=w[j];j++)
for(uint t=pri[i];1LL*t*pri[i]<=w[j];t*=pri[i])
sum[j]+=sum[get(w[j]/t)]+f[i]*(g[get(w[j]/t)]-i+1);
}
/*
inline uint S(uint x,int y){
if(x<=1 || pri[y]>x) return 0;
uint ret=0;
for(int i=y;i<=cnt && 1LL*pri[i]*pri[i]<=x;i++){
uint p1=pri[i],p2=1LL*pri[i]*pri[i];
for(uint t=pri[i];t*pri[i]<=x;t*=pri[i])
ret+=S(x/t,i+1)+f[i]*(g[get(x/t)]-(i-1));
} return ret;
}
*/
inline uint solve(uint x){
int pos=get(x);
if(h[pos]!=-1) return h[pos];
uint ret=sum[pos]+g[pos];//minp>=pri[1] & prime
for(uint i=2,j;i<=x;i=j+1){
j=x/(x/i); ret-=solve(x/i)*(j-i+1);
} return h[pos]=ret;
}
int main(){
rd(n); rd(K); sqr=sqrt(n);
memset(h,-1,sizeof h);
get_prime(sqr); init();
uint lst=0,cur,x;
for(uint i=1,j;i<=n;i=j+1){
x=n/i,j=n/x; cur=solve(j);
ans+=(cur-lst)*x*x; lst=cur;
}
printf("%u\n",ans);
return 0;
}
想一想
m
i
n
_
25
min\_25
min_25筛就是一些套路的东西，只要能推出式子就什么都好办了
Q
A
Q
QAQ
QAQ关键是怎么推式子啊