【jzoj1938】【2011集训队出题】【Crash的数字表格】【莫比乌斯反演】

题目大意

Ans=∑Nx=1∑My=1Lcm(x,y)

解题思路

Ans=∑Nx=1∑My=1Lcm(x,y)

=∑Nx=1∑My=1xy/Gcd(x,y)

=∑Nd=11/d∑Nx=1∑My=1xy[Gcd(x,y)==d] 默认N<=M

g[d]=∑Nx=1∑My=1xy[Gcd(x,y)==d]

f[d]=∑Nx=1∑My=1xy[d|Gcd(x,y)]=d2(1+⌊N/d⌋)⌊N/d⌋(1+⌊M/d⌋)⌊M/d⌋/4=∑⌊N/d⌋i=1g[di]

根据莫比乌斯反演

g[d]=∑⌊N/d⌋i=1f[di]mu(i)

=∑⌊N/d⌋i=1d2i2(1+⌊N/di⌋)⌊N/di⌋(1+⌊M/di⌋)⌊M/di⌋/4mu(i)

Ans=∑Nd=11/d∑⌊N/d⌋i=1d2i2(1+⌊N/di⌋)⌊N/di⌋(1+⌊M/di⌋)⌊M/di⌋/4mu(i)

令T=id

=∑NT=1T(1+⌊N/T⌋)⌊N/T⌋(1+⌊M/T⌋)⌊M/T⌋/4∑i|Tmu(i)i

由莫比乌斯函数性质得 h[i]=∑i|Tmu(i)i 是积性函数， h[x]=∑pi1−pi 可以用线筛求。

对于前一部分可以分块求，对于后一部分可以前缀和求。

前面的T其实可以放进后面求。

code

#include<cmath> #include<cstdio> #include<cstring> #include<algorithm> #define LF double #define LL long long #define Min(a,b) ((a<b)?a:b) #define Max(a,b) ((a>b)?a:b) #define Fo(i,j,k) for(int i=j;i<=k;i++) #define Fd(i,j,k) for(int i=j;i>=k;i--) using namespace std; int const Mxn=1e7,Mo=20101009; int N,M,F[Mxn+9],Pri[Mxn+9]; int main(){ freopen("d.in","r",stdin); freopen("d.out","w ",stdout); scanf("%d%d",&N,&M); if(N>M)swap(N,M); F[1]=1; Fo(i,2,N){ if(!F[i])Pri[++Pri[0]]=i,F[i]=1-i; Fo(j,1,Pri[0]){ if(i*Pri[j]>N)break; if(i%Pri[j]==0){F[i*Pri[j]]=F[i];break;} else F[i*Pri[j]]=F[i]*(1-Pri[j]); } F[i]=(1ll*F[i]*i+F[i-1])%Mo; } int Ans=0; for(int i=1,j,Ni,Mi;i<=N;i=j+1){ Ni=N/i;Mi=M/i; j=Min(N/(Ni),M/(Mi)); Ans=(Ans+((1ll+Ni)*Ni/2%Mo)*((1ll+Mi)*Mi/2%Mo)%Mo %Mo*(F[j]-F[i-1])%Mo)%Mo; } printf("%d",(Ans+Mo)%Mo); return 0; }