Min_25筛学习笔记
0、一些定义
- 定义域为正整数集，陪域为复数域的函数被称为数论函数。
- 设\(f(x)\)是一个数论函数，若对于任意一对互质的正整数\(a,b\)，有\(f(ab)=f(a)f(b)\)，则称\(f(x)\)是积性函数。
- 记\(\delta(x)\)为\(x\)的最小质因子，\(\gamma(x)\)为\(x\)的最大质因子。
- 记\(cnt(x)\)为\(x\)的质因子个数。
- 记第\(i\)个质数为\(p_i\)。特别地，\(p_0=1\)。
- 记\(\pi(x)\)为不超过\(x\)的质数个数。由素数定理，\(\pi(x)=O(\frac{x}{\log x})\)。
1、min_25筛
论文题：定义\(\sigma_0(n)\)为\(n\)的正约数个数，求\(S(n)=\sum_{i=1}^n\sigma_0(i^k)\)。\(n,k\leq 10^{10}\)。
拓展：给定一个积性函数\(f(x)\)，且对于任意质数\(p\)和正整数\(e\)，满足：
- \(f(p)\)是一个关于\(p\)的低次多项式。
- 计算\(f(p^e)\)的复杂度可忽略。
- 求\(S(n)=\sum_{i=1}^nf(i)\)。
解法
令\(f(x)=\sigma_0(x^k)\)，容易发现\(f(x)\)是积性函数。
令
那么\(S(n)=f(1)+g(n,0)\)。
考虑计算\(g(x)\)的递推式，将\(x\)分成质数和合数两个部分，合数部分枚举最小质因子(一定\(\leq \sqrt{n}\))，质数部分单独计算，于是
令
那么
如果能快速计算\(h\)，就可以在很低的复杂度里递推求出\(g(n,0)\)。注意到\(f\)在质数处的取值为多项式，所以只需要对于非负整数\(s\)，计算\(h_s(n)=\sum_{1\leq i\leq \pi(n)}p_i^s\)。注意到\(g(n,i)\)的递推过程，我们只需要对于\(1,2,\ldots,\lfloor\sqrt{n}\rfloor,\left\lfloor\frac{n}{\lfloor\sqrt{n}\rfloor}\right\rfloor,\ldots,\left\lfloor\frac{n}{2}\right\rfloor,n\)这\(O(\sqrt{n})\)个数求出\(h\)值。
令\(h_s'(n,i)\)，为前\(n\)个数，经过\(i\)次埃氏筛后剩下的数的\(s\)次幂之和。那么\(h_s(n)=h_s'(n,\pi(\sqrt{n}))\)。
对于计算\(h_s'(n,i)\)的递推式，考虑第\(i\)轮埃氏筛的过程：
- 若\(n<p_i^2\)，那么不会有任何数被筛去，于是\(h_s'(n,i)=h_s'(n,i-1)\)。
- 若\(n\geq p_i^2\)，那么被筛去的数必有质因子\(p_i\)，于是应当减去\(p_i^sh_s'\left(\left\lfloor\frac{n}{p_i}\right\rfloor,i-1\right)\)。但是会有最小质因子小于\(p_i\)的数被减去，所以应当加上这部分的贡献，为\(p_i^sh_s'(p_{i-1},i-1)\)。
那么
边界条件为
暴力实现即可。
时间复杂度
\(O(\frac{n^\frac{3}{4}}{\log n})\)
证明咕着。
例题
注意\(f\)在质数\(p\)上的取值为\(p^2-p\)，于是需要同时求出质数的和和平方和，即\(h_1,h_2\)。然后照着公式实现即可。
#include<cstdio>
#define ll long long
#define P 1000000007
#define inv3 333333336
#define N 200005
ll n;
int prm[N],_prm;
bool vis[N];
inline void sieve(int x){
for(int i=2;i<=x;i++){
if(!vis[i])
prm[++_prm]=i;
for(int j=1;j<=_prm&&i*prm[j]<=x;j++){
vis[i*prm[j]]=1;
if(!(i%prm[j]))
break;
}
}
}
int sq,m;
inline int id(ll x){
return x<=sq?x:m-n/x+1;
}
ll num[N];
int h1[N],h2[N];//h(x)=h2(x)-h1(x)，h1,h2开成滚动数组
inline void get_h(){
while(1ll*sq*sq<=n)
sq++;
sq--;
sieve(sq);
for(int i=1;i<=sq;i++)
num[++m]=i;
for(int i=sq;i;i--)
if(n/i>sq)
num[++m]=n/i;
for(int i=1;i<=m;i++){
int tmp=num[i]%P;
h1[i]=1ll*tmp*(tmp+1)/2%P-1;
h2[i]=1ll*tmp*(tmp+1)/2%P*(2*tmp+1)%P*inv3%P-1;
}//边界条件
for(int i=1;i<=_prm;i++)
for(int j=m;num[j]>=1ll*prm[i]*prm[i];j--){//注意倒序递推
int k=id(num[j]/prm[i]);
h1[j]=(h1[j]-1ll*prm[i]*(h1[k]-h1[prm[i-1]]+P)%P+P)%P;
h2[j]=(h2[j]-1ll*prm[i]*prm[i]%P*(h2[k]-h2[prm[i-1]]+P)%P+P)%P;
}
}
inline int g(ll x,int i){
int k=id(x),res=((h2[k]-h1[k]+P)%P-(h2[prm[i]]-h1[prm[i]]+P)%P+P)%P;//h(x)-h(p_i)
for(int j=i+1;j<=_prm&&1ll*prm[j]*prm[j]<=x;j++)
for(ll pe=prm[j];pe*prm[j]<=x;pe=pe*prm[j]){
int tmp=pe%P;
res=(res+1ll*tmp*(tmp-1)%P*g(x/pe,j)%P)%P;//f(p_j^e)*g(x/p_j^e,j)
res=(res+1ll*tmp*prm[j]%P*(1ll*tmp*prm[j]%P-1)%P)%P;//f(p_j^(e+1))
}
return res;
}
int main(){
scanf("%lld",&n);
get_h();
printf("%d\n",((g(n,0))+1)%P);
}
