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

题目大意

$Ans=\sum_{x=1}^{N}\sum_{y=1}^{M}Lcm(x,y)$

解题思路

$Ans=\sum_{x=1}^{N}\sum_{y=1}^{M}Lcm(x,y)$

$=\sum_{x=1}^{N}\sum_{y=1}^{M}xy/Gcd(x,y)$

$=\sum_{d=1}^{N}1/d\sum_{x=1}^{N}\sum_{y=1}^{M}xy[Gcd(x,y)==d]$默认N<=M

$g[d]=\sum_{x=1}^{N}\sum_{y=1}^{M}xy[Gcd(x,y)==d]$

$f[d]=\sum_{x=1}^{N}\sum_{y=1}^{M}xy[d|Gcd(x,y)]=d^2(1+\lfloor N/d\rfloor)\lfloor N/d\rfloor (1+\lfloor M/d\rfloor)\lfloor M/d\rfloor /4=\sum_{i=1}^{\lfloor N/d\rfloor}g[di]$

根据莫比乌斯反演

$g[d]=\sum_{i=1}^{\lfloor N/d\rfloor}f[di]mu(i)$

$=\sum_{i=1}^{\lfloor N/d\rfloor}d^2i^2(1+\lfloor N/di\rfloor)\lfloor N/di\rfloor (1+\lfloor M/di\rfloor)\lfloor M/di\rfloor /4mu(i)$

$Ans=\sum_{d=1}^N1/d\sum_{i=1}^{\lfloor N/d\rfloor}d^2i^2(1+\lfloor N/di\rfloor)\lfloor N/di\rfloor (1+\lfloor M/di\rfloor)\lfloor M/di\rfloor /4mu(i)$

令T=id

$=\sum_{T=1}^NT(1+\lfloor N/T\rfloor)\lfloor N/T\rfloor (1+\lfloor M/T\rfloor)\lfloor M/T\rfloor /4\sum_{i|T}mu(i)i$

由莫比乌斯函数性质得$h[i]=\sum_{i|T}mu(i)i$是积性函数，$h[x]=\sum_{p_i}1-p_i$可以用线筛求。

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

前面的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;
}