关于min_25筛的一些简单笔记

90nwyn

2020-10-07 13:32:32

Personal

设$f(i)$是一个积性函数,满足$f(p),p \in prime$可以用多项式表示,同时$f(p^k)$可以快速计算出来 考虑如何求解$\sum _{i=1}^nf(i)(1<n<10^{11})$ ------------ 另一个新的问题:对于积性函数$f(i)$,求$\sum_{i=1}^n[i \in prime] f(i)$ 为了方便理解,问题变为求$\sum_{i=1}^n[i \in prime] i^k$ 定义$prime_i$为第$i$个质数,$minp(i)$为$i$的最小质因子,令$minp(1)=inf$ 对于$\sqrt n$范围内的前$i$个质数的$k$次方前缀和$sum_i=\sum _{j=1}^i prime_j^k$可以线性筛解决 设$g(n,i)=\sum_{j=1}^n[j \in prime $ $or$ $ minp(j)>prime_i]j^k$ 直观地来看,$g(n,i)$就是在 $n$ 范围内的埃式筛算法进行第$i$轮后尚未被筛去的数的 $k$ 次方和 显然,$g(n,0)=\sum_{i=2}^ni^k$ 思考从$g(n,i-1)$到$g(n,i)$的转移 - $prime_i^2>n$,在埃氏筛的第 $i$ 轮没有数字被筛掉,$g(n,i)=g(n,i-1)$ - $prime_i^2 \leq n$,考虑$prime_i$在$g(n,i-1)$中的贡献,$g(n,i)=g(n,i-1)-prime_i^k(g(\left \lfloor \frac{n}{prime_i}\right \rfloor,i-1)-sum_{i-1})$ 令$cnt$为$\sqrt n$范围内的质数个数,有 $\sum_{i=1}^n[i \in prime] i^k=sum_{cnt}+g(n,cnt)$ 时间复杂度 $O(\frac{n^{3/4}}{logn})$ ------------ 现在问题回归到求解$\sum _{i=1}^nf(i)$ 定义$s(n,i)=\sum_{j=1}^n[minp(j) \geq prime_i] f(j)$ 那么,$\sum _{i=1}^nf(i)=s(n,1)+f(1)$ 问题最终化为如何求解$s(n,i)$ 因为$\sum_{i=1}^n[i \in prime] f(i)$已经可以快速计算,那么不难得到答案中质数的贡献$\sum_{j=1}^n[j \in prime] f(j)-\sum_{j=1}^{i-1}f(prime_j)$ 考虑合数的贡献,枚举每个合数的最小质因子以及最小质因子的次幂,由$f$为积性函数得 $\sum_{j=1}^{prime_j^2 \leq n}\sum_{k=1}^{prime_j^{k+1} \leq n} ( s( \left \lfloor \frac{n}{prime_j^k}\right \rfloor,j+1)f(prime_j^k)+f(prime_j^{k+1}))$ ------------ [题目链接](https://vjudge.net/problem/LibreOJ-6053) 考虑$f(p),p \in prime$的值 若$p=2$,$f(p)=p-1$,反之$f(p)=p+1$ 问题相当于min_25筛求出$n$以内的质数和与质数个数 ------------ ```cpp #include<bits/stdc++.h> using namespace std; typedef long long int ll; const int M=1e6+5,mod=1e9+7,inv2=500000004; int prime[M],tot,vis[M],sqr,sum[M],id1[M],id2[M],m,g[M],h[M]; ll n,w[M]; void init(int siz) { for(int i=2;i<=siz;i++) { if(!vis[i])prime[++tot]=i,sum[tot]=(sum[tot-1]+i)%mod; for(int j=1;j<=tot&&i*prime[j]<=siz;j++) { vis[i*prime[j]]=1; if(i%prime[j]==0)break; } } } int getid(ll x){return (x<=sqr)?id1[x]:id2[n/x];} ll s(ll x,int y) { if(x<=1||prime[y]>x)return 0; int pos=getid(x); ll res=(g[pos]-sum[y-1]-h[pos]+y-1+mod)%mod; if(y==1)res+=2; for(int i=y;i<=tot&&(ll)prime[i]*prime[i]<=x;i++) { ll p1=prime[i]; ll p2=p1*p1; for(int e=1;p2<=x;e++,p1=p2,p2*=prime[i]) res=(res+s(x/p1,i+1)*(prime[i]^e)%mod+(prime[i]^(e+1))%mod)%mod; } return res; } int main() { scanf("%lld",&n); sqr=sqrt(n); init(sqr); for(ll l=1,r=0;l<=n;l=r+1) { w[++m]=n/l; r=n/w[m]; if(w[m]<=sqr)id1[w[m]]=m; else id2[n/w[m]]=m; g[m]=((2+w[m])%mod)*((w[m]-1)%mod)%mod*inv2%mod; h[m]=(w[m]-1)%mod; } for(int j=1;j<=tot;j++) for(int i=1;i<=m&&(ll)prime[j]*prime[j]<=w[i];i++) { int k=getid(w[i]/prime[j]); g[i]=(g[i]-(ll)prime[j]*((g[k]-sum[j-1]+mod)%mod)+mod)%mod; h[i]=(h[i]-h[k]+j-1+mod)%mod; } printf("%lld\n",(s(n,1)+1)%mod); return 0; } ```