本文深入探讨了Min_25筛算法，一种高效计算积性函数前缀和的方法。通过定义和推导，文章详细解释了如何利用最小素因子筛选技术优化计算过程，并给出了具体的实现代码。适用于算法竞赛和高级数学应用。
关于min_25筛：
现在有一个积性函数，对素数 p p p有 f ( p ) = a p i f(p)=ap^i f(p)=api，我想求出他的前缀和。
设
m
i
n
p
(
x
)
minp(x)
minp(x)表示
x
x
x最小的素因子。
设
q
(
n
,
j
)
=
∑
i
=
1
n
[
m
i
n
p
(
x
)
≥
p
j
]
f
(
i
)
q(n,j)=\sum_{i=1}^n[minp(x)\ge p_j]f(i)
q(n,j)=∑i=1n[minp(x)≥pj]f(i)，即最小素因子不小于第
j
j
j个素数
p
j
p_j
pj的
f
f
f值的和。
由于
q
(
n
,
j
)
q(n,j)
q(n,j)比
q
(
n
,
j
+
1
)
q(n,j+1)
q(n,j+1)多出来的只有
m
i
n
p
(
x
)
=
p
j
minp(x)=p_j
minp(x)=pj的部分，有：
q
(
n
,
j
)
=
q
(
n
,
j
+
1
)
+
∑
i
=
1
p
j
i
≤
n
q
(
⌊
n
p
j
i
⌋
,
j
+
1
)
f
(
p
j
i
)
+
∑
i
=
1
p
j
i
≤
n
f
(
p
j
i
)
q(n,j)=q(n,j+1)+\sum_{i=1}^{p_j^i\le n}q(\lfloor \frac n {p_j^i}\rfloor,j+1)f(p_j^i)+\sum_{i=1}^{p_j^i\le n}f(p_j^i)
q(n,j)=q(n,j+1)+i=1∑pji≤nq(⌊pjin⌋,j+1)f(pji)+i=1∑pji≤nf(pji)
注意到
q
(
n
,
j
)
q(n,j)
q(n,j)中，最小的合数为
p
j
2
p_j^2
pj2，故对于
p
j
2
>
n
p_j^2>n
pj2>n的情况，有：
q
(
n
,
j
)
=
∑
i
=
j
p
i
≤
n
f
(
p
i
)
q(n,j)=\sum_{i=j}^{p_i\le n}f(p_i)
q(n,j)=i=j∑pi≤nf(pi)
设
g
(
n
,
j
)
=
∑
i
=
1
n
[
i
∈
p
r
i
m
e
∨
m
i
n
p
(
i
)
>
p
j
]
f
(
i
)
g(n,j)=\sum_{i=1}^n[i\in prime \lor minp(i)>p_j]f(i)
g(n,j)=∑i=1n[i∈prime∨minp(i)>pj]f(i)，即不大于
n
n
n的所有素数的
f
f
f值和最小素因子大于
p
j
p_j
pj的合数的
f
f
f值的和。
设
z
z
z为最大的满足
p
z
2
≤
n
p_z^2\le n
pz2≤n的值，即第
z
z
z个素数是最大的不大于
⌊
n
⌋
\lfloor \sqrt n\rfloor
⌊n⌋的素数，则
g
(
n
,
p
z
)
g(n,pz)
g(n,pz)表示不大于
n
n
n的所有素数的
f
f
f值的和。
由于
g
(
n
,
j
)
g(n,j)
g(n,j)比
g
(
n
,
j
−
1
)
g(n,j-1)
g(n,j−1)减少的只有
m
i
n
p
(
x
)
=
p
j
minp(x)=p_j
minp(x)=pj的合数的部分，有：
g
(
n
,
j
)
=
g
(
n
,
j
−
1
)
−
(
g
(
⌊
n
p
j
⌋
,
j
−
1
)
−
∑
i
=
1
j
−
1
f
(
p
i
)
)
f
(
p
j
)
g(n,j)=g(n,j-1)-(g(\lfloor \frac n {p_j} \rfloor,j-1)-\sum_{i=1}^{j-1}f(p_i))f(p_j)
g(n,j)=g(n,j−1)−(g(⌊pjn⌋,j−1)−i=1∑j−1f(pi))f(pj)
同时，当
p
j
2
>
n
p_j^2>n
pj2>n时，容易得：
g
(
n
,
j
)
=
g
(
n
,
j
−
1
)
g(n,j)=g(n,j-1)
g(n,j)=g(n,j−1)
以此，我们便可以计算出
p
j
2
>
n
p_j^2>n
pj2>n时的
q
q
q值，即：
q
(
n
,
j
)
=
g
(
n
,
j
)
−
∑
i
=
1
j
−
1
f
(
p
i
)
=
g
(
n
,
p
z
)
−
∑
i
=
1
j
−
1
f
(
p
i
)
∣
n
≥
p
j
q
(
n
,
j
)
=
0
∣
n
<
p
j
q(n,j)=g(n,j)-\sum_{i=1}^{j-1}f(p_i)=g(n,pz)-\sum_{i=1}^{j-1}f(p_i)\ |n\ge p_j\\ q(n,j)=0\ |n < p_j
q(n,j)=g(n,j)−i=1∑j−1f(pi)=g(n,pz)−i=1∑j−1f(pi) ∣n≥pjq(n,j)=0 ∣n<pj
（需要注意的是，在
g
g
g的递推式中我们将
f
f
f当作一个完全积性函数，此时我们设
g
g
g的边界
g
(
n
,
0
)
=
∑
i
=
1
n
f
(
i
)
g(n,0)=\sum_{i=1}^nf(i)
g(n,0)=∑i=1nf(i)，容易证明最终计算出的
g
(
n
,
p
z
)
g(n,pz)
g(n,pz)是正确的。故当
f
f
f是一个包含常数个项的多项式时，需要将每一项分别计算，否则结果将是错误的，除非能够证明
f
f
f本身就是一个完全积性函数。）
计算
g
g
g的复杂度为：（我tm怎么知道怎么算)
∑
i
=
1
n
n
l
n
n
+
n
i
l
n
n
i
=
O
(
n
3
4
l
n
n
)
\sum_{i=1}^{\sqrt n}\frac {\sqrt n} {ln {\sqrt n}}+\frac {\sqrt {\frac n i}} {ln{\sqrt {\frac n i}}}=O(\frac {n^{\frac 3 4}} {lnn})
i=1∑nlnnn+lninin=O(lnnn43)
计算
q
q
q的复杂度为：（因为连式子都不会写所以直接写答案好了）
O
(
n
3
4
l
n
n
)
O(\frac {n^{\frac 3 4}} {lnn})
O(lnnn43)
题目
枚举gcd，有：
a
n
s
=
∑
i
=
1
n
(
i
m
i
n
p
(
i
)
)
k
F
(
⌊
n
i
⌋
)
F
(
x
)
=
∑
i
=
1
x
∑
j
=
1
x
[
(
i
,
j
)
=
1
]
1
ans=\sum_{i=1}^n(\frac i {minp(i)})^kF(\lfloor \frac n i\rfloor)\\ F(x)=\sum_{i=1}^x\sum_{j=1}^x[(i,j)=1]1
ans=i=1∑n(minp(i)i)kF(⌊in⌋)F(x)=i=1∑xj=1∑x[(i,j)=1]1
F
F
F很好求，考虑如何求
(
i
m
i
n
p
(
i
)
)
k
(\frac i {minp(i)})^k
(minp(i)i)k的前缀和。
设：
g
(
n
)
=
∑
i
=
1
n
[
i
∈
p
r
i
m
e
]
i
g
′
(
n
)
=
∑
i
=
1
n
[
i
∈
p
r
i
m
e
]
i
k
q
(
n
,
j
)
=
∑
i
=
1
n
[
m
i
n
p
(
i
)
≥
p
j
]
i
k
q
′
(
n
,
j
)
=
∑
i
=
1
n
[
m
i
n
p
(
i
)
≥
p
j
]
(
i
m
i
n
p
(
i
)
)
k
g(n)=\sum_{i=1}^n[i\in prime]i\\ g'(n)=\sum_{i=1}^n[i\in prime]i^k\\ q(n,j)=\sum_{i=1}^n[minp(i)\ge p_j]i^k\\ q'(n,j)=\sum_{i=1}^n[minp(i)\ge p_j](\frac i {minp(i)})^k
g(n)=i=1∑n[i∈prime]ig′(n)=i=1∑n[i∈prime]ikq(n,j)=i=1∑n[minp(i)≥pj]ikq′(n,j)=i=1∑n[minp(i)≥pj](minp(i)i)k
有：
q
(
n
,
j
)
=
q
(
n
,
j
+
1
)
+
∑
d
=
1
p
j
d
≤
n
(
1
+
q
(
⌊
n
p
j
d
⌋
,
j
+
1
)
)
p
j
d
k
q
′
(
n
,
j
)
=
q
′
(
n
,
j
+
1
)
+
∑
d
=
1
p
j
d
≤
n
(
1
+
q
(
⌊
n
p
j
d
⌋
,
j
+
1
)
)
p
j
(
d
−
1
)
k
q
′
(
n
,
j
)
=
M
A
X
{
0
,
g
(
n
)
−
j
+
1
}
∣
p
j
2
>
n
q
(
n
,
j
)
=
M
A
X
{
0
,
g
′
(
n
)
−
∑
i
=
1
j
−
1
p
i
}
∣
p
j
2
>
n
q(n,j)=q(n,j+1)+\sum_{d=1}^{p_j^d\le n}(1+q(\lfloor \frac n {p_j^d}\rfloor,j+1))p_j^{dk}\\ q'(n,j)=q'(n,j+1)+\sum_{d=1}^{p_j^d\le n}(1+q(\lfloor \frac n {p_j^d}\rfloor,j+1))p_j^{(d-1)k}\\ q'(n,j)=MAX\{0,g(n)-j+1\}|p_j^2>n\\ q(n,j)=MAX\{0,g'(n)-\sum_{i=1}^{j-1}p_i\}|p_j^2>n
q(n,j)=q(n,j+1)+d=1∑pjd≤n(1+q(⌊pjdn⌋,j+1))pjdkq′(n,j)=q′(n,j+1)+d=1∑pjd≤n(1+q(⌊pjdn⌋,j+1))pj(d−1)kq′(n,j)=MAX{0,g(n)−j+1}∣pj2>nq(n,j)=MAX{0,g′(n)−i=1∑j−1pi}∣pj2>n
g
g
g和
g
′
g'
g′可以用min_25筛计算，
q
q
q和
q
′
q'
q′直接递推即可（当然他们长得其实和min_25筛一样）。
至于如何在不能用逆元的情况下求
∑
i
=
1
n
i
k
\sum_{i=1}^ni^k
∑i=1nik，也就是
g
′
(
n
,
0
)
g'(n,0)
g′(n,0)：
设
s
(
n
,
k
)
s(n,k)
s(n,k)表示将
n
n
n个不同的球放进
k
k
k个相同的箱子的方案数，即第二类斯特林数。
由其组合性质，有：
n
k
=
∑
i
=
1
n
(
n
i
)
s
(
k
,
i
)
i
!
∑
i
=
1
n
i
k
=
∑
i
=
1
n
∑
j
=
1
i
(
i
j
)
s
(
k
,
j
)
j
!
=
∑
j
=
1
n
s
(
k
,
j
)
j
!
∑
k
=
1
n
(
k
j
)
=
∑
j
=
1
n
s
(
k
,
j
)
j
!
(
k
+
1
j
+
1
)
=
∑
j
=
1
n
s
(
k
,
j
)
(
k
+
1
)
j
+
1
‾
j
+
1
n^k=\sum_{i=1}^n\binom n i s(k,i)i!\\ \sum_{i=1}^ni^k=\sum_{i=1}^n\sum_{j=1}^i\binom i js(k,j)j!\\ =\sum_{j=1}^ns(k,j)j!\sum_{k=1}^n\binom k j\\ =\sum_{j=1}^ns(k,j)j!\binom {k+1} {j+1}\\ =\sum_{j=1}^ns(k,j)\frac {(k+1)^{\underline {j+1}}} {j+1}
nk=i=1∑n(in)s(k,i)i!i=1∑nik=i=1∑nj=1∑i(ji)s(k,j)j!=j=1∑ns(k,j)j!k=1∑n(jk)=j=1∑ns(k,j)j!(j+1k+1)=j=1∑ns(k,j)j+1(k+1)j+1
这样就可以无视逆元
O
(
k
2
)
O(k^2)
O(k2)求出结果了。
#include<algorithm>
#include<iostream>
#include<cstring>
#include<climits>
#include<cstdio>
#include<cmath>
using namespace std;
#define clr(a) memset(a,0,sizeof(a))
//--Container
//--
typedef unsigned int ui;typedef long long ll;
const int up=1e6;const int dup=5e4;
ui strl[51][51];
bool bd[up+10];int n,dk,pn,pr[up+10],moi[up+10];ui prs[dup+10],moip[up+10];
void _init(){
int i,j,k,d,t;for(pn=0,moi[1]=1,i=2;i<=up;++i){
if(!bd[i]){pr[++pn]=i;moi[i]=i-1;}
for(j=1;j<=pn&&(ll)pr[j]*i<=up;++j){
bd[pr[j]*i]=1,moi[pr[j]*i]=moi[i]*(pr[j]-1);
if(!(i%pr[j])){
moi[pr[j]*i]=moi[i]*pr[j];break;
}
}
}
for(strl[0][0]=1,i=1;i<=50;++i)for(j=1;j<=i;++j)
strl[i][j]=strl[i-1][j-1]+strl[i-1][j]*j;
for(moip[1]=1,i=2;i<=up;++i)moip[i]+=moip[i-1]+(moi[i]<<1);
};
ui _qni(ui a,int b){
ui r=1;for(;b;b>>=1,a*=a)if(b&1)r*=a;
return r;
};
ui _dom(ui n,ui i){
ui rs=1;int j;for(j=0;j<i;++j,--n)rs*=((n%i)?n:n/i);
return rs;
};
ui _prnm(int n){
ui rs=0;int i,j;for(i=1;i<=min(dk,n);++i)rs+=strl[dk][i]*_dom(n+1,i+1);
return rs;
};
ui psm[(dup<<1)+10],pssm[(dup<<1)+10];int id[2][dup+10],mg[(dup<<1)+10],dn,pz,dt;
inline int _id(int x){return x<=dn?id[0][x]:id[1][n/x];};
void _genid(){
int b,e,z;for(dn=sqrt(n),dt=0,b=1;b<=n;b=e+1){
e=n/(n/b);z=n/e;if(z<=dn)id[0][z]=++dt;
else
id[1][n/z]=++dt;
mg[dt]=z;
}
};
void _ps(){
int i,j,k,d,t,z;for(i=1;i<=pn&&pr[i]*pr[i]<=n;++i)prs[i]+=prs[i-1]+_qni(pr[i],dk);pz=i-1;
for(i=1;i<=dt;++i)psm[i]=mg[i]-1,pssm[i]=_prnm(mg[i])-1;
for(i=1;i<=pz;++i)for(z=pr[i]*pr[i],j=1;j<=dt&&z<=mg[j];++j){
psm[j]-=(psm[_id(mg[j]/pr[i])]-i+1);
pssm[j]-=(prs[i]-prs[i-1])*(pssm[_id(mg[j]/pr[i])]-prs[i-1]);
}
};
ui dsm[(dup<<1)+10],dssm[(dup<<1)+10];
inline ui _dsm(int x,int j){
return pr[j]*pr[j]>x?(x>pr[j-1]?pssm[_id(x)]-prs[j-1]:0):dsm[_id(x)];
};
inline ui _dssm(int x,int j){
return pr[j]*pr[j]>x?(x>pr[j-1]?psm[_id(x)]-j+1:0):dssm[_id(x)];
};
void _ds(){
int i,j,z;ui t,d;ll k;for(i=pz;i;--i)for(z=pr[i]*pr[i],j=1;j<=dt&&z<=mg[j];++j){
t=prs[i]-prs[i-1],k=pr[i],d=1;dsm[j]=_dsm(mg[j],i+1),dssm[j]=_dssm(mg[j],i+1);
for(;k<=mg[j];t*=(prs[i]-prs[i-1]),k*=pr[i],d*=(prs[i]-prs[i-1])){
dsm[j]+=_dsm(mg[j]/k,i+1)*t+t;
dssm[j]+=_dsm(mg[j]/k,i+1)*d+d;
}
}
};
ui mem[dup+10];
inline ui _mem(int x){return x<=up?moip[x]:mem[_id(x)];};
void __mem(int x){
int i,j,b,e,z,t;ui rs=(ui)x*x;for(b=2;b<=x;b=e+1){
e=x/(x/b);rs-=_mem(x/e)*(e-b+1);
}
mem[_id(x)]=rs;
};
void _me(){
int i,j,k,d,t;for(i=1;mg[i]>up;++i);for(--i;i;__mem(mg[i--]));
};
void cl(){
int i,j,k,b,e;scanf("%d %d",&n,&dk);ui rs=0;
_init();_genid();_ps();_ds();_me();
for(b=2;b<=n;b=e+1){
e=n/(n/b);rs+=(_dssm(e,1)-_dssm(b-1,1))*_mem(n/e);
}
printf("%u\n",rs);
};
int main(){
#ifndef ONLINE_JUDGE
freopen("in.txt","r",stdin);
freopen("out.txt","w",stdout);
#endif
cl();
return 0;
};
扫一扫
348
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
