标签：25 const min int 质数 St mod sum 小记
min_25筛
由 dalao min_25 发明的筛子，据说时间复杂度是极其优秀的 $ O(\frac {n^{\frac 3 4}} {\log n}) $，常数还小。
1. 质数 $ k $ 次方前缀和（基础）
求 $ \sum_{p \leq n}p^k $
我们考虑一个 $ \rm DP $ 的思路：设 $ g(n,j) $ 为：
\[\sum_{i=1}^n[(\sum_{t=1}^j[p_t|i])=0] i^k \]其实就是不大于 $ n $ 的，且不含有 $ p_1 $ ~ $ p_j $ 中任意一个质因子的数的 $ k $ 次方之和。
我们考虑从 $ g(n,j-1) $ 转移上来。
发现 $ g(n,j) $ 相对于 $ g(n,j-1) $ 减少的就是含有 $ p_k $ 质因子的数的 $ k $ 次方之和，于是我们得到了转移方程：
\[g(n,k) = g(n,k-1)-p_j^kg(\lfloor \frac n {p_j} \rfloor,j) \]但是，这个方程似乎有点儿问题？
我们把 $ p_j $ 以内的质数算重了，所以要加上。
最终的方程就是：
\[g(n,k) = g(n,k-1)-p_j^k(g(\lfloor \frac n {p_j} \rfloor,j)-g(p_{j-1},j-1) \]不过如果我们要维护 $ O(n) $ 个 $ g $ 复杂度开销是很大的，所以考虑 玄学优化 去掉不重要的东西。
首先我们注意到，在 $ p_j > \sqrt n $ 是，对答案的贡献只有 $ p_j^k $，所以我们可以先只考虑 $ p_j \leq \sqrt n $ 的情况。
$ g(p_{j-1},j-1) $ 显然可以在 $ O(\sqrt n) $ 的复杂度内线性筛预处理出来，就只需要考虑 $ g(\lfloor \frac n {p_j} \rfloor,j) $
然后根据整除分块的结论，我们只需要处理 $ O(\sqrt n) $ 个 $ g $ 就行了。。。
然后复杂度就降下来了，我也不会证，不过看上去似乎是 $ O(\sqrt n\log\log n) $ 的（？）
但是 $ n $ 太大需要离散化。。。
所以我们依照杜教筛的离散化，把大于 $ \sqrt n $ 的映射到 $ 1 $ ~ $ \sqrt n $。
给一下参考代码（这里只算了质数的 $ 1 $ 次方和 $ 2 $ 次方之和）：
const int St=1e6+5,mod=1e9+7;
inline int Add(const int&a,const int&b){
return a+b>=mod?a+b-mod:a+b;
}
inline int Get(const ll&k){
return k<=sqr?id1[k]:id2[n/k];
}
inline void seive(const int&n){
int i,j,x;zhi[1]=1;
for(i=2;i<=n;++i){
if(!zhi[i]){
pri[++top]=i;
sum1[top]=Add(sum1[top-1],i);
sum2[top]=Add(sum2[top-1],1ll*i*i%mod);
}
for(j=1;j<=top&&(x=i*pri[j])<=n;++j){
zhi[x]=1;
if(!(i%pri[j]))break;
}
}
}
ll min_25(const ll&n){
int i,j;
seive(sqr=sqrt(n));
for(ll L=1,R;L<=n;L=R+1){
R=n/(n/L);
w[++tot]=n/L;
g1[tot]=w[tot]%mod;
g2[tot]=g1[tot]*(g1[tot]+1)/2%mod*(g1[tot]*2+1)%mod*inv3%mod-1;
g1[tot]=g1[tot]*(g1[tot]+1)/2%mod-1;
if(n/L<=sqr)id1[n/L]=tot;
else id2[n/(n/L)]=tot;
}
for(i=1;i<=top;++i){
for(j=1;j<=tot&&pri[i]*pri[i]<=w[j];++j){
ll k=Get(w[j]/pri[i]);
g1[j]-=pri[i]*(g1[k]-sum1[i-1]+mod)%mod;
g2[j]-=pri[i]*pri[i]%mod*(g2[k]-sum2[i-1]+mod)%mod;
if(g1[j]<=0)g1[j]+=mod;
if(g2[j]<=0)g2[j]+=mod;
}
}
}
在实现的时候将 $ g $ 的第二维滚掉了qwq
2.min_25筛
没错前面的都是铺垫（雾）
让我们看看我们要求的问题：
\[\sum_{i=1}^nf(i) \]其中 $ f $ 是一个低阶多项式，其在质数乘方处的取值能够快速算出。
将答案分成两部分考虑：
\[\sum_{p \leq n}f(p)+\sum_{p^e \leq n,p \leq \sqrt n}f(p^e)(\sum_{i=1 \And minp > p}^{\lfloor \frac n {p^e} \rfloor}f(i)) \]其中 $ minp $ 是 $ i $ 最小的质因子。
按照类似上面 $ \rm DP $ 的套路，设 $ S(n,x) $ 是不大于 $ n $ 的，且不含有 $ p_1 $ ~ $ p_x $ 中任何一个因子的 $ f $ 之和。
显然答案是 $ S(n,0) $。
再将答案分成两部分，一部分是大于 $ p_x $ 的质数的 $ f $ 之和，另一部分就是剩下的。
第一部分，也就是 $ g(n) - g(p_x) $。
对照着上面的式子，我们可以得到以下的转移方程
\[S(n,x) = g(n) - g(p_x) +\sum_{p_k^e \leq n,k > x}f(p_k^e)(S(\lfloor \frac n {p_k^e} \rfloor,k)+[e!=1]) \]然后还有一些细节，就是关于 $ 1 $ 的问题。。。
由于是照着 $ \rm wucstdio $ dalao 的题解学的，所以我也没有算 $ 1 $，而是在结束时加上。
例题
\[p^k(p^k-1) = p^{2k} - p^k \]在质数处就是 $ p^2 - p $
只需要筛出质数的 $ 1 $ 次方前缀和和 $ 2 $ 次方前缀和即可。
不要问我为什么上面的代码筛的也是 $ 1 $ 次和 $ 2 $ 次前缀和，我懒总行了吧
#include<cstdio>
#include<cmath>
typedef long long ll;
const int St=1e6+5,mod=1e9+7,inv3=333333336;
int sqr,top,tot,zhi[St],id1[St],id2[St];
ll n,w[St],g1[St],g2[St],pri[St],sum1[St],sum2[St];
inline int Add(const int&a,const int&b){
return a+b>=mod?a+b-mod:a+b;
}
inline int Get(const ll&k){
return k<=sqr?id1[k]:id2[n/k];
}
inline void seive(const int&n){
int i,j,x;zhi[1]=1;
for(i=2;i<=n;++i){
if(!zhi[i]){
pri[++top]=i;
sum1[top]=Add(sum1[top-1],i);
sum2[top]=Add(sum2[top-1],1ll*i*i%mod);
}
for(j=1;j<=top&&(x=i*pri[j])<=n;++j){
zhi[x]=1;
if(!(i%pri[j]))break;
}
}
}
ll S(const ll&n,const int&k){
if(pri[k]>=n)return 0;
ll x=Get(n),ans=Add((g2[x]-g1[x]+sum1[k]-sum2[k])%mod,mod);
for(int i=k+1;i<=top&&pri[i]*pri[i]<=n;++i){
ll p=pri[i];
for(int e=1;p<=n;++e,p*=pri[i]){
ll id=p%mod;
ans=Add(ans,id*(id-1)%mod*((e!=1)+S(n/p,i))%mod);
}
}
return ans;
}
ll min_25(const ll&n){
int i,j;
seive(sqr=sqrt(n));
for(ll L=1,R;L<=n;L=R+1){
R=n/(n/L);
w[++tot]=n/L;
g1[tot]=w[tot]%mod;
g2[tot]=g1[tot]*(g1[tot]+1)/2%mod*(g1[tot]*2+1)%mod*inv3%mod-1;
g1[tot]=g1[tot]*(g1[tot]+1)/2%mod-1;
if(n/L<=sqr)id1[n/L]=tot;
else id2[n/(n/L)]=tot;
}
for(i=1;i<=top;++i){
for(j=1;j<=tot&&pri[i]*pri[i]<=w[j];++j){
ll k=Get(w[j]/pri[i]);
g1[j]-=pri[i]*(g1[k]-sum1[i-1]+mod)%mod;
g2[j]-=pri[i]*pri[i]%mod*(g2[k]-sum2[i-1]+mod)%mod;
if(g1[j]<=0)g1[j]+=mod;
if(g2[j]<=0)g2[j]+=mod;
}
}
return Add(S(n,0),1);
}
signed main(){
scanf("%lld",&n);
printf("%d",min_25(n));
}
然后是一道很经典的三倍经验题DIVCNTK，做出这道题后 DIVCNT2 和 DIVCNT3 就是三倍经验。
设积性函数 $ f(p^e) = d((pe)k) = d(p^{ek}) $
在质数处的取值就是 $ k+1 $
#include<cstdio>
#include<cmath>
typedef unsigned long long ull;
const int M=2e6+5;
ull n,m,lim,sqr,tot,top,w[M],g[M],pri[M],zhi[M],id1[M],id2[M];
inline int Get(const ull&k){
return k<=sqr?id1[k]:id2[n/k];
}
inline void sieve(const int&n){
int i,j,x;zhi[1]=1;
for(i=2;i<=n;++i){
if(!zhi[i])pri[++top]=i;
for(j=1;j<=top&&(x=i*pri[j])<=n;++j){
zhi[x]=1;
if(!(i%pri[j]))break;
}
}
}
ull S(const ull&n,const int&k){
if(pri[k]>=n)return 0;
ull x=Get(n),ans=g[x]-k*(m+1);
for(ull i=k+1;i<=lim&&pri[i]*pri[i]<=n;++i){
for(ull e=1,p=pri[i];p<=n;++e,p*=pri[i]){
ans+=(e*m+1)*(S(n/p,i)+(e!=1));
}
}
return ans;
}
ull min_25(const ull&n){
int i,j;
sqr=sqrt(n);tot=0;lim=1;
while(pri[lim]*lim[pri]<=n)++lim;--lim;
for(ull L=1,R;L<=n;L=R+1){
R=n/(n/L);
w[++tot]=n/L;
g[tot]=w[tot]-1;
if(n/L<=sqr)id1[n/L]=tot;
else id2[n/(n/L)]=tot;
}
for(i=1;i<=lim;++i){
for(j=1;j<=tot&&pri[i]*pri[i]<=w[j];++j){
ull k=Get(w[j]/pri[i]);
g[j]-=g[k]-i+1;
}
}
for(i=1;i<=tot;++i)g[i]*=(m+1);
return S(n,0)+1;
}
signed main(){
int i,T;
sieve(2e6);
scanf("%d",&T);
for(i=1;i<=T;++i){
scanf("%llu%llu",&n,&m);
printf("%llu\n",min_25(n));
}
}
标签：25,const,min,int,质数,St,mod,sum,小记 来源： https://www.cnblogs.com/lmpp/p/15721158.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。