P1829 [国家集训队]Crash的数字表格 / JZPTAB

传送门

$ Ans = sum_i^Nsum_j^Mlcm(i,j)$

因为 $lcm(x,y)=xy/gcd(x,y)$

所以 $Ans = sum_i^Nsum_j^Mij/gcd(i,j)=sum_i^Nsum_j^M( i/gcd(i,j) )( j/gcd(i,j) ) cdot gcd(i,j)$

改成枚举所有可能的 $gcd(i,j)=t$，变成 $Ans=sum_tsum_i^{left lfloor frac{N}{t} ight floor}sum_j^{left lfloor frac{M}{t} ight floor}ijt[gcd(i,j)==1]$

显然莫比乌斯反演

设 $f[x]=sum_tsum_i^{left lfloor frac{N}{t} ight floor}sum_j^{left lfloor frac{M}{t} ight floor}ijt[gcd(i,j)==x]$，$F[x]=sum_tsum_i^{left lfloor frac{N}{t} ight floor}sum_j^{left lfloor frac{M}{t} ight floor}ijt[x|gcd(i,j)]$

那么 $Ans=f[1]$

那么有 $F[x]=sum_{x|d}f[d]$，直接反演

    $f[1]=sum_{1|d}mu(d/1)sum_tsum_i^{left lfloor frac{N}{t} ight floor}sum_j^{left lfloor frac{M}{t} ight floor}ijt[d|gcd(i,j)]$

$=sum_dmu(d)sum_tsum_i^{left lfloor frac{N}{td} ight floor}sum_j^{left lfloor frac{M}{td} ight floor}ijtd^2$

枚举 $T=td$，则原式 $=sum_dmu(d)sum_{d|T}sum_i^{left lfloor frac{N}{T} ight floor}sum_j^{left lfloor frac{M}{T} ight floor}ijTd$

$=(sum_dmu(d)d)(sum_{d|T}T)(sum_i^{left lfloor frac{N}{T} ight floor}i)(sum_j^{left lfloor frac{M}{T} ight floor}j)$

改成先枚举 $T$ 再枚举 $d$

$=(sum_TT)(sum_{d|T}mu(d)d)(sum_i^{left lfloor frac{N}{T} ight floor}i)(sum_j^{left lfloor frac{M}{T} ight floor}j)$

显然后面的 $(sum_i^{left lfloor frac{N}{T} ight floor}i)$ 和 $(sum_j^{left lfloor frac{M}{T} ight floor}j)$ 可以预处理

考虑能不能预处理 $g[T]=sum_{d|T}mu(d)d$

发现好像可以线性筛，考虑假设此时已经知道了 $g[x]$，且 $y=xp$，$p$ 为 $y$ 唯一一个最小的质因数

那么 $g[xp]=sum_{d|xp}mu(d)d=(sum_{d|x}mu(d)d)+(sum_{d|x}mu(dp)dp)$

因为 $p$ 为 $y$ 唯一一个最小的质因数，所以 $gcd(d,p)=1$ ，所以 $mu(dp)=-mu(d)$

则 $g[xp]=(sum_{d|x}mu(d)d)-(sum_{d|x}mu(d)dp)=(sum_{d|x}mu(d)d)(1-p)=(1-p)g[x]$

如果 $y=xp$，$p$ 为 $y$ 最小的一个质因数但不唯一，则根据 $g$ 的定义发现 $p$ 对 $x$ 没有贡献，所以此时 $g[xp]=g[x]$

然后直接线性筛就行了，具体看代码

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
using namespace std;
typedef long long ll;
inline int read()
{
    int x=0,f=1; char ch=getchar();
    while(ch<'0'||ch>'9') { if(ch=='-') f=-1; ch=getchar(); }
    while(ch>='0'&&ch<='9') { x=(x<<1)+(x<<3)+(ch^48); ch=getchar(); }
    return x*f;
}
const int N=2e7+7,mo=20101009;
inline int fk(int x) { return x>=mo ? x-mo : x; }
int n,m,mi,pri[N],mu[N],g[N],tot;
bool not_pri[N];
ll ans;
void pre()
{
    mu[1]=1; g[1]=1;
    for(int i=2;i<=mi;i++)
    {
        if(!not_pri[i]) pri[++tot]=i,mu[i]=mo-1,g[i]=mo+1-i;
        for(int j=1;j<=tot;j++)
        {
            ll t=1ll*i*pri[j]; if(t>mi) break;
            not_pri[t]=1; if(!(i%pri[j])) { g[t]=g[i]; break; }
            mu[t]=mo-mu[i]; g[t]= 1ll*g[i]*(mo+1-pri[j])%mo;
        }
    }
}
inline int sum(int x) { return 1ll*x*(x+1)/2%mo; }
int main()
{
    n=read(),m=read();
    mi=min(n,m); pre();
    for(int _d=1;_d<=mi;_d++)
    {
        ll t=g[_d],s=1ll*sum(n/_d)*sum(m/_d)%mo;
        ans=fk(ans+ ((t*s%mo)*_d%mo) );
    }
    printf("%lld
",(ans%mo+mo)%mo);
    return 0;
}