51nod 1056

时间:2023-03-09 20:36:19
51nod 1056

n<=10000000000

51nod 1056

然后欧拉函数的前缀和可以用莫比乌斯函数的前缀和快速求,注意各种取模

#include<cstdio>
typedef long long i64;
const int N=,P=,M=,I2=(M+)/;
int ps[N/],pp=;
bool isnp[N+];
int ms[N+],phs[N+];
i64 xs[P];
int ys[P];
int f(i64 n){
if(n<=N)return ms[n];
int w=n%P;
while(xs[w]){
if(xs[w]==n)return ys[w];
w+=;
if(w>=P)w-=P;
}
xs[w]=n;
i64 s=;
for(i64 i=,j;i<=n;i=j+){
j=n/(n/i);
s+=(j-i+)*f(n/i);
}
return ys[w]=-s;
}
int sum_phi(i64 n){
if(n<=N)return phs[n];
i64 ans=;
for(i64 i=,j;i<=n;i=j+){
i64 k=n/i;
j=n/k;
k%=M;
ans=(ans+k*k%M*(f(j)-f(i-)))%M;
}
ans=(ans+)*I2%M;
return ans;
}
int main(){
ms[]=phs[]=;
for(int i=;i<=N;i++){
if(!isnp[i])ms[ps[pp++]=i]=-,phs[i]=i-;
for(int j=;j<pp&&i*ps[j]<=N;j++){
isnp[i*ps[j]]=;
if(i%ps[j])ms[i*ps[j]]=-ms[i],phs[i*ps[j]]=phs[i]*(ps[j]-);
else{
phs[i*ps[j]]=phs[i]*ps[j];
break;
}
}
}
for(int i=;i<=N;i++)ms[i]+=ms[i-];
for(int i=;i<=N;i++)(phs[i]+=phs[i-])%=M;
i64 n,ans=;
scanf("%lld",&n);
for(i64 i=,j=,k;i<=n;i=j+){
k=n/i;
j=n/k;
k%=M;
ans=(ans+k*(k+)%M*(sum_phi(j)-sum_phi(i-)))%M;
}
n%=M;
ans=(ans-n*(n+)%M*I2)%M;
printf("%lld",(ans+M)%M);
return ;
}