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

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

题目

解法

首先简化题意

\[ans=\sum_{i=1}^N\sum_{j=1}^Mlcm(i,j) \]

\[ans=\sum_{i=1}^N\sum_{j=1}^M\frac{ij}{gcd(i,j)} \]

枚举一下\(gcd\)

\[ans=\sum_{k=1}^N\frac{1}{k}\sum_{i=1}^N\sum_{j=1}^M{ij}[gcd(i,j)=k] \]

把\(k\)提出来

\[ans=\sum_{k=1}^Nk\sum_{i=1}^\frac{N}{k}\sum_{j=1}^\frac{M}{k}{ij}[gcd(i,j)=1] \]

后面那一部分莫比乌斯反演走起
定义

\[f(x)=\sum_{k=1}^Nk\sum_{i=1}^\frac{N}{k}\sum_{j=1}^\frac{M}{k}{ij}[gcd(i,j)=x] \]

\[g(x)=\sum_{x|d}f(d) \]

所以\(ans=f(1)\)

\[g(x)=\sum_{x|d}\sum_{k=1}^Nk\sum_{i=1}^\frac{N}{k}\sum_{j=1}^\frac{M}{k}{ij}[gcd(i,j)=d] \]

去掉第一个\(\sum\)

\[g(x)=\sum_{k=1}^Nk\sum_{i=1}^\frac{N}{k}\sum_{j=1}^\frac{M}{k}{ij}[gcd(i,j)\%x=0] \]

提出\(x\)，排除\(gcd\)的影响

\[g(x)=\sum_{k=1}^Nkx^2\sum_{i=1}^\frac{N}{kx}\sum_{j=1}^\frac{M}{kx}{ij} \]

根据莫比乌斯定理

\[f(x)=\sum_{x|d}μ(\frac{d}{x})g(d) \]

\[f(x)=\sum_{x|d}μ(\frac{d}{x})\sum_{k=1}^Nkd^2\sum_{i=1}^\frac{N}{kd}\sum_{j=1}^\frac{M}{kd}{ij} \]

而\(ans=f(1)\)

\[ans=\sum_{d=1}^Nμ(d)\sum_{k=1}^Nkd^2\sum_{i=1}^\frac{N}{kd}\sum_{j=1}^\frac{M}{kd}ij \]

\[ans=\sum_{d=1}^Nμ(d)\sum_{k=1}^Nkd^2\sum_{i=1}^\frac{N}{kd}i\sum_{j=1}^\frac{M}{kd}j \]

后面两个是等差数列求和

\[ans=\sum_{d=1}^Nμ(d)\sum_{k=1}^Nkd^2\frac{(\frac{N}{kd})(\frac{N}{kd}+1)(\frac{M}{kd})(\frac{M}{kd}+1)}{4} \]

(我一开始做的时候是在这里就结束了，两次数论分块时间复杂度变\(O(n)\),但后来做了一道一样的但是有多组询问的题目，所以就有了以下的优化)

设\(D=kd\)，枚举\(D\),考虑到若\(D>N\)，对答案没有贡献,所以\(D\)只要枚举到\(N\)

\[ans=\sum_{D=1}^ND\sum_{d|D}μ(d)d\frac{(\frac{N}{D})(\frac{N}{D}+1)(\frac{M}{D})(\frac{M}{D}+1)}{4} \]

中间那个线性筛一下
然后前缀和加数论分块是\(O(\sqrt n)\)的。

完整代码(未优化的)(下面有优化过的)

#include<bits/stdc++.h>
#define ll long long
#define il inline
#define rg register
using namespace std;
const ll N=1e7+20,p=20101009;
__int128 ljq=1;
ll pre[N],prime[N],n,m,s;
bool isprime[N];
void work(){
    cin>>n>>m;if(n>m)swap(n,m);
    for(rg ll i=2;i<=n;++i){
        if(!isprime[i])pre[i]=-1,prime[++s]=i;
        for(rg ll j=1;j<=s&&i*prime[j]<=n;++j){
            isprime[i*prime[j]]=1;
            if(i%prime[j])pre[i*prime[j]]=-pre[i];
            else break;
        }
    }
    pre[1]=1;
    for(int i=1;i<=n;++i)
        pre[i]=(ljq*pre[i]*i*i+pre[i-1])%p;
}
ll js(rg int x,rg int y){
    ll ans=0;
    for(rg int l=1,r;l<=x;l=r+1){
        r=min(x/(x/l),y/(y/l));
        ans+=(ljq*(pre[r]-pre[l-1])*(x/l+1)*(x/l)/2*(y/l)*(y/l+1)/2)%p;
        ans%=p;
    }
    return ans;
}
int main(){
    work();
    ll ans=0;
    for(int l=1,r;l<=n;l=r+1){
        r=min(n/(n/l),m/(m/l));
        ans+=(ljq*(l+r)*(r-l+1)/2*js(n/l,m/l))%p;
        ans%=p;
    }
    cout<<((ans%p)+p)%p<<endl;
    return 0;
}

优化代码(注意模数不一样)

#include<bits/stdc++.h>
#define ll long long
#define rg register
using namespace std;
const int N=1e7+20,yl=1e8+9,M=1e7+10;
ll dis[N],isprime[N],cnt,prime[N];
void eular(){
	for(rg int i=2;i<=M;++i){
		if(!isprime[i])prime[++cnt]=i,dis[i]=1-i;
		for(int j=1;j<=cnt&&i*prime[j]<=M;++j){
			isprime[i*prime[j]]=1;
			if(i%prime[j])dis[i*prime[j]]=dis[i]*dis[prime[j]]%yl;
			else {dis[i*prime[j]]=dis[i];break;}
		}
	}
	dis[1]=1;
	for(int i=1;i<=M;++i){
		dis[i]=(dis[i]*i+dis[i-1])%yl;
	}
}
int main(){
	eular();
	int T;
	cin>>T;
	while(T--){
		int n,m;
		scanf("%d%d",&n,&m);
		if(n>m)swap(n,m);
		ll ans=0;
		for(ll l=1,r;l<=n;l=r+1){
			r=min(n/(n/l),m/(m/l));
			ans=(ans+(dis[r]-dis[l-1])*((n/l)*(n/l+1)/2%yl)%yl*((m/l)*(m/l+1)/2%yl))%yl;
		}
		if(ans<0)ans+=yl;
		printf("%lld\n",ans);
	}
}