51Nod 1584 加权约数和
51Nod 1584 加權(quán)約數(shù)和
要求:\(\sum_{i=1}^n\sum_{j=1}^n\max(i,j)\sigma_1(ij)\)
枚舉\(\max(i,j)\),式子變?yōu)?#xff1a;\(2\sum_{i=1}^ni\sum_{j=1}^i\sigma_1(ij)-\sum_{i=1}^ni\sigma_1(i^2)\)
這里不加證明地給出結(jié)論:\(\sigma_1^k(ij)=\sum_{p|i}\sum_{q|j}(p\frac{j}{q})^k[\gcd(p,q)=1]\)
先看右邊部分,可以線性篩\(\sigma_1(n^2)\)。
(篩約數(shù)個(gè)數(shù)和的方法:\(\sigma_1(n=\prod_{i=1}^mp_i^{e_i})=\prod_{i=1}^m(1+p_i+\cdots+p_i^{e_i})\),篩平方的約數(shù)個(gè)數(shù)和類似,事實(shí)上任意\(k\)次方的都能篩)
或者也可以推式子:
\(\sum_{i=1}^ni\sigma_1(i^2)\)
\(=\sum_{i=1}^ni\sum_{p|i}\sum_{q|i}p\frac iq\sum_{d|p,d|q}\mu(d)\)
\(=\sum_{i=1}^ni\sum_{p|i}\sum_{q|i}p\frac iq\sum_{d|p,d|q}\mu(d)\)
\(=\sum_{i=1}^ni\sum_{d|i}\mu(d)\sum_{d|p}^ip\sum_{d|q}^i\frac iq\)
\(=\sum_{i=1}^ni\sum_{d|i}\mu(d)\sum_{p|\frac id}^idp\sum_{q|\frac id}^i\frac i{dq}\)
\(=\sum_{i=1}^ni\sum_{d|i}\mu(d)d\sigma_1(\frac id)^2\)
左邊部分類似,直接寫(xiě)結(jié)論了:
\(\sum_{i=1}^ni\sum_{j=1}^i\sigma_1(ij)=\sum_{i=1}i\sum_{d|i}\mu(\frac id)(\frac id)\sigma_1(d)S_{\sigma_1}(d)\)
其中\(S_f(n)=\sum_{i=1}^nf(i)\)
右邊用線性篩:
#include<bits/stdc++.h> #define il inline #define vd void #define mod 1000000007 typedef long long ll; il ll gi(){ll x=0,f=1;char ch=getchar();while(!isdigit(ch))f^=ch=='-',ch=getchar();while(isdigit(ch))x=x*10+ch-'0',ch=getchar();return f?x:-x; } int pr[1000010],P,sd[1000010],_sd[1000010],sum[1000010],ssd[1000010],mu[1000010]; int sd2[1000010],_sd2[1000010],ssd2[1000010],sum2[1000010]; bool yes[1000010]; int ans[1000010]; int main(){ #ifdef XZZSBfreopen("in.in","r",stdin);freopen("out.out","w",stdout); #endifint N=1000000;sd[1]=sd2[1]=1;mu[1]=1;for(int i=2;i<=N;++i){if(!yes[i])pr[++P]=i,sd[i]=sum[i]=i+1,_sd[i]=1,mu[i]=mod-1,sd2[i]=sum2[i]=(1+i+1ll*i*i)%mod,_sd2[i]=1;for(int j=1;j<=P&&i*pr[j]<=N;++j){yes[i*pr[j]]=1;if(i%pr[j]==0){sum[i*pr[j]]=(1ll*sum[i]*pr[j]+1)%mod;sd[i*pr[j]]=1ll*_sd[i]*sum[i*pr[j]]%mod;_sd[i*pr[j]]=_sd[i];sum2[i*pr[j]]=(1ll*sum2[i]*pr[j]%mod*pr[j]+pr[j]+1)%mod;sd2[i*pr[j]]=1ll*_sd2[i]*sum2[i*pr[j]]%mod;_sd2[i*pr[j]]=_sd2[i];mu[i*pr[j]]=0;break;}sum[i*pr[j]]=1+pr[j];sd[i*pr[j]]=1ll*sd[i]*sum[i*pr[j]]%mod;_sd[i*pr[j]]=sd[i];sum2[i*pr[j]]=(1ll*pr[j]*pr[j]+pr[j]+1)%mod;sd2[i*pr[j]]=1ll*sd2[i]*sum2[i*pr[j]]%mod;_sd2[i*pr[j]]=sd2[i];mu[i*pr[j]]=mod-mu[i];}}for(int i=1;i<=N;++i)ssd[i]=(ssd[i-1]+sd[i])%mod;for(int i=1;i<=N;++i)ssd2[i]=(ssd2[i-1]+1ll*i*sd2[i])%mod;for(int i=1;i<=N;++i)for(int j=i;j<=N;j+=i)if(mu[i])ans[j]=(ans[j]+2ll*j*mu[i]%mod*i%mod*sd[j/i]%mod*ssd[j/i])%mod;for(int i=1;i<=N;++i)ans[i]=(ans[i]+ans[i-1])%mod;int T=gi(),n;for(int i=1;i<=T;++i)n=gi(),printf("Case #%d: %d\n",i,(ans[n]-ssd2[n]+mod)%mod);return 0; }右邊用推式子:
#include<bits/stdc++.h> #define il inline #define vd void #define mod 1000000007 typedef long long ll; il ll gi(){ll x=0,f=1;char ch=getchar();while(!isdigit(ch))f^=ch=='-',ch=getchar();while(isdigit(ch))x=x*10+ch-'0',ch=getchar();return f?x:-x; } int pr[1000010],P,sd[1000010],_sd[1000010],sum[1000010],ssd[1000010],mu[1000010]; bool yes[1000010]; int ans[1000010]; int main(){ #ifdef XZZSBfreopen("in.in","r",stdin);freopen("out.out","w",stdout); #endifint N=1000000;sd[1]=1;mu[1]=1;for(int i=2;i<=N;++i){if(!yes[i])pr[++P]=i,sd[i]=i+1,_sd[i]=1,sum[i]=i+1,mu[i]=mod-1;for(int j=1;j<=P&&i*pr[j]<=N;++j){yes[i*pr[j]]=1;if(i%pr[j]==0){sum[i*pr[j]]=(1ll*sum[i]*pr[j]+1)%mod;sd[i*pr[j]]=1ll*_sd[i]*sum[i*pr[j]]%mod;_sd[i*pr[j]]=_sd[i];mu[i*pr[j]]=0;break;}sum[i*pr[j]]=1+pr[j];sd[i*pr[j]]=1ll*sd[i]*sum[i*pr[j]]%mod;_sd[i*pr[j]]=sd[i];mu[i*pr[j]]=mod-mu[i];}}for(int i=1;i<=N;++i)ssd[i]=(ssd[i-1]+sd[i])%mod;for(int i=1;i<=N;++i)for(int j=i;j<=N;j+=i){ans[j]=(ans[j]+2ll*j*mu[j/i]%mod*(j/i)%mod*sd[i]%mod*ssd[i])%mod;ans[j]=(ans[j]-1ll*j*mu[i]%mod*i%mod*sd[j/i]%mod*sd[j/i]%mod+mod)%mod;}for(int i=1;i<=N;++i)ans[i]=(ans[i]+ans[i-1])%mod;int T=gi(),n;for(int i=1;i<=T;++i)n=gi(),printf("Case #%d: %d\n",i,ans[n]);return 0; }轉(zhuǎn)載于:https://www.cnblogs.com/xzz_233/p/11142149.html
總結(jié)
以上是生活随笔為你收集整理的51Nod 1584 加权约数和的全部?jī)?nèi)容,希望文章能夠幫你解決所遇到的問(wèn)題。
- 上一篇: Linux网络性能评估工具iperf 、
- 下一篇: ajax轮播图控件,基于json数据的j