SCUT - 114 - 作业之数学篇 - 杜教筛
\(A(n)=\sum\limits_{i=1}^{n} \frac{lcm(i,n)}{gcd(i,n)}\)
\(=\sum\limits_{i=1}^{n} \frac{in}{gcd^2(i,n)}\)
枚举g:
\(A(n)=n\sum\limits_{g|n}\frac{1}{g^2} \sum\limits_{i=1}^{n} i [gcd(i,n)==g]\)
最内层除以g:
\(A(n)=n\sum\limits_{g|n}\frac{1}{g} \sum\limits_{i=1}^{\frac{n}{g}} i [gcd(i,\frac{n}{g})==1]\)
考虑里面的:
\(B(n)=\sum\limits_{i=1}^{n} i [gcd(i,n)==1]\)
显然:
\(B(n)=\frac{1}{2}n([n==1]+\varphi(n))\)
代回 \(A(n)\) 里面:
\(A(n)=\sum\limits_{g|n}\frac{n}{g} B(\frac{n}{g})\)
\(=\sum\limits_{g|n} g B(g)\)
\(=\frac{1}{2}(1+\sum\limits_{g|n} g^2 \varphi(g))\)
代回 \(S(n)\) 里面:
$S(n)=\sum\limits_{i=1}^{n} A(i) \(
\)= \sum\limits_{i=1}^{n} \frac{1}{2}(1+\sum\limits_{g|i} g^2 \varphi(g))\(
\)= \frac{n}{2}+\frac{1}{2}\sum\limits_{i=1}^{n}\sum\limits_{g|i} g^2 \varphi(g)$
记:
\(C(n)=\sum\limits_{i=1}^{n}\sum\limits_{g|i} g^2 \varphi(g)\)
显然:
\(C(n)=\sum\limits_{d=1}^{n}d^2 \varphi(d) \lfloor\frac{n}{d}\rfloor\)
后面是一个分块,前面是个杜教筛。然后我们召唤鹏哥就可以通过这道题了。
转化为快速求 $\sum\limits_{i=1}^{n} i^2 \varphi(i) $,卷一个id平方然后查鹏哥的cheatsheet过了。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod=1e9+7;
const int inv2=(mod+1)>>1;
const int MAXN=4e6;
int pri[MAXN+1];
int &pritop=pri[0];
ll B[MAXN+1];
ll PrefixB[MAXN+1];
int pk[MAXN+1];
void sieve(int n=MAXN) {
pk[1]=1;
B[1]=1;
for(int i=2; i<=n; i++) {
if(!pri[i]) {
pri[++pritop]=i;
pk[i]=i;
ll tmp=1ll*i*i;
//.
if(tmp>=mod)
tmp%=mod;
B[i]=tmp*(i-1);
if(B[i]>=mod)
B[i]%=mod;
//.
}
for(int j=1; j<=pritop; j++) {
int &p=pri[j];
int t=i*p;
//.
if(t>n)
break;
pri[t]=1;
if(i%p) {
pk[t]=p;
B[t]=B[i]*B[p];
if(B[t]>=mod)
B[t]%=mod;
//.
} else {
pk[t]=pk[i]*p;
if(pk[t]==t) {
ll tmp=1ll*t*t;
if(tmp>=mod)
tmp%=mod;
B[t]=tmp*(t-i);
if(B[t]>=mod)
B[t]%=mod;
//.
} else {
B[t]=B[t/pk[t]]*B[pk[t]];
if(B[t]>=mod)
B[t]%=mod;
}
break;
}
}
}
for(int i=1; i<=n; i++) {
PrefixB[i]=PrefixB[i-1]+B[i];
if(PrefixB[i]>=mod)
PrefixB[i]-=mod;
}
}
inline int phi(int n) {
int res=n;
for(int i=1; i<=pritop&&pri[i]*pri[i]<=n; i++) {
if(n%pri[i]==0) {
res-=res/pri[i];
while(n%pri[i]==0)
n/=pri[i];
}
if(n==1)
return res;
}
res-=res/n;
return res;
//.
}
inline ll qpow(ll x,int n) {
ll res=1;
while(n) {
if(n&1) {
res*=x;
if(res>=mod)
res%=mod;
}
x*=x;
if(x>=mod)
x%=mod;
n>>=1;
}
return res;
}
int inv4=qpow(4,mod-2);
int inv6=qpow(6,mod-2);
inline ll s2(ll n) {
ll tmp=n*(n+1);
if(tmp>=mod)
tmp%=mod;
tmp*=(n*2+1);
if(tmp>=mod)
tmp%=mod;
tmp*=inv6;
if(tmp>=mod)
tmp%=mod;
return tmp;
}
inline ll s3(ll n) {
ll tmp=n*(n+1);
if(tmp>=mod)
tmp%=mod;
tmp*=tmp;
if(tmp>=mod)
tmp%=mod;
tmp*=inv4;
if(tmp>=mod)
tmp%=mod;
return tmp;
}
unordered_map<int,ll> PB;
inline ll S(int n) {
if(n<=MAXN)
return PrefixB[n];
if(PB.count(n))
return PB[n];
ll ret=s3(n);
for(int l=2,r; l<=n; l=r+1) {
int t=n/l;
r=n/t;
ll tmp=s2(r)-s2(l-1);
if(tmp<0)
tmp+=mod;
tmp*=S(t);
if(tmp>=mod)
tmp%=mod;
ret-=tmp;
if(ret<0)
ret+=mod;
}
return PB[n]=ret;
}
inline ll Prefix(int l,int r) {
l--;
ll PL=S(l);
ll PR=S(r);
ll res=PR-PL;
if(res<0)
res+=mod;
return res;
}
inline ll P(int n) {
ll res=0;
for(int l=1,r; l<=n; l=r+1) {
int t=n/l;
r=n/t;
ll tmp=1ll*t*Prefix(l,r);
if(tmp>=mod)
tmp%=mod;
res+=tmp;
if(res>=mod)
res-=mod;
}
return res;
}
inline ll Ans(int n) {
ll res=(1ll*n*inv2)%mod;
res+=(P(n)*inv2)%mod;
return res%mod;
}
int main() {
#ifdef Yinku
freopen("Yinku.in","r",stdin);
#endif // Yinku
sieve();
int n;
while(~scanf("%d",&n)) {
printf("%lld\n",Ans(n));
}
return 0;
}