[整理]min_25 筛
大概是一种牛逼的筛法,本质思想是 dp。
方便起见先规定本文的一些符号:
\(p_i\) 表示第 \(i\) 个质数,\(p\) 表示任意一个质数;
\(\text{P}\) 表示全体质数;
\(\min_i\) 表示 \(i\) 的最小质因数。
现在我们要求一个积性函数 \(f\) 的前缀和,这个 \(f\) 需要满足在 \(p^k\) 处为关于 \(p\) 的多项式,这样我们只考虑 \(f(p)=p^k\) 的情况(是一个完全积性函数)。
容易发现如果设 \(S(n,j)\) 为最小质因数大于 \(p_j\) 的 \(f\) 之和,我们要求的就是 \(S(n,0)\),而下面可以看到这个 \(S\) 是更加好求的。
接下来就是比较牛逼的了,我们设一个 dp 数组 \(g_{n,j}\) 表示 \(\sum\limits_{i=1}^n[i\in\text{P}\lor\min_i>p_j]i^k\)。
显然 \(p_{j-1}^2>n\) 时这个东西就等价于 \(g_{n,j-1}\),下面考虑其他情况。
容易发现比起 \(g_{n,j-1}\) 少的部分就是最小质因数恰好为 \(p_j\) 的合数。利用完全积性函数的性质,可以提一个 \(p_j^k\) 出来,剩下的就是 \(g_{\lfloor n/p_j\rfloor,j-1}\)。
吗?我们多减了一些质数,应该加回来一个 \(g_{p_{j-1},j-1}\)。
求出了 \(g\),我们尝试寻找 \(g\) 和 \(S\) 之间的关系。
显然,\(S(n,j)\) 可以按照类似的处理方式提取质因数后递归地计算:\(S(n,j)=g_{n,\text{P}}-\sum\limits_{i=1}^{j-1}f(p_i)+\sum\limits_{i>j,p_i^x\le n}f(p_i^x)(S(\lfloor n/p_i^x\rfloor,i)+[x>1])\)。
然后因为变量都是 \(n\) 整除一个数,可以离散化后做到空间 \(\sqrt{n}\) 级别。
非常丑的洛谷模板题代码:
const int N=400010,p=1e9+7;
int n,B,g1[N],g2[N],id1[N],id2[N],inv3;
bool np[N];int pri[N],sp1[N],sp2[N],pri0;
il int Pow(int a,int b=p-2){
int res=1;
for(;b;a=(LL)a*a%p,b>>=1)if(b&1)res=(LL)res*a%p;
return res;
}
il void Init(){
np[1]=1;
for(int i=2;i<N;i++){
if(!np[i]){
pri[++pri0]=i;
sp1[pri0]=(sp1[pri0-1]+i)%p;
sp2[pri0]=(sp2[pri0-1]+i*i%p)%p;
}
for(int j=1;j<=pri0&&i*pri[j]<N;j++){
np[i*pri[j]]=1;
if(i%pri[j]==0)break;
}
}
}
int w[N],tot;
il void DP(){
//j=0 时就是求所有数的和
for(int l=1,r;l<=n;l=r+1){
r=n/(n/l);
w[++tot]=n/l,g1[tot]=w[tot]%p;
g2[tot]=g1[tot]*(g1[tot]+1)/2%p*(2*g1[tot]+1)%p*inv3%p;
g1[tot]=(g1[tot]*(g1[tot]+1)/2%p-1+p)%p;
g2[tot]=(g2[tot]-1+p)%p;
if(n/l<=B)id1[n/l]=tot;
else id2[n/(n/l)]=tot;
}
for(int i=1;i<=pri0;i++){
for(int j=1;j<=tot&&pri[i]*pri[i]<=w[j];j++){
int qwq=w[j]/pri[i];
if(qwq<=B)qwq=id1[qwq];
else qwq=id2[n/qwq];
g1[j]=(g1[j]-pri[i]*(g1[qwq]-sp1[i-1]+p)%p+p)%p;
g2[j]=(g2[j]-pri[i]*pri[i]%p*\
(g2[qwq]-sp2[i-1]+p)%p+p)%p;
}
}
}
int S(int x,int j){
if(pri[j]>=x)return 0;
int qwq=x<=B?id1[x]:id2[n/x];
int res=(g2[qwq]-g1[qwq]-(sp2[j]-sp1[j]+p)%p+p+p)%p;
for(int i=j+1;i<=pri0&&pri[i]*pri[i]<=x;i++){
int pw=pri[i];
for(int k=1;pw<=x;k++,pw=pw*pri[i]){
int fuck=pw%p;
(res+=fuck*(fuck-1)%p*(S(x/pw,i)+(k>1))%p)%=p;
}
}
return res;
}
signed main(){
Read(n),inv3=Pow(3),B=sqrt(n);
Init(),DP();
printf("%lld\n",(S(n,0)+1)%p);
KafuuChino HotoKokoa
}