题意:
给出n,mn,mn,m,计算
∑i=1n∑j=1mlcm(i,j)
\sum_{i=1}^{n}\sum_{j=1}^{m} lcm(i,j)
i=1∑nj=1∑mlcm(i,j)
Solution:
lcmlcmlcm不好处理转化为gcdgcdgcd,有
∑i=1n∑j=1mijgcd(i,j)
\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{ij}{gcd(i,j)}
i=1∑nj=1∑mgcd(i,j)ij
优先枚举因子gcdgcdgcd,这样可以转化成整除问题
∑d=1n∑i=1n∑j=1m[gcd(i,j)=d]ijd
\sum_{d=1}^{n}\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]\frac{ij}{d}
d=1∑ni=1∑nj=1∑m[gcd(i,j)=d]dij
利用gcd(i,j)=d⇒gcd(id,jd)=1gcd(i,j)=d\Rightarrow gcd(\frac{i}{d},\frac{j}{d})=1gcd(i,j)=d⇒gcd(di,dj)=1,枚举右式的id,jd\frac{i}{d},\frac{j}{d}di,dj作为新的i,ji,ji,j
∑d=1n∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]id⋅jdd=∑d=1n∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]ijd
\sum_{d=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]\frac{id\cdot jd}{d}=\sum_{d=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]ijd
d=1∑ni=1∑⌊dn⌋j=1∑⌊dm⌋[gcd(i,j)=1]did⋅jd=d=1∑ni=1∑⌊dn⌋j=1∑⌊dm⌋[gcd(i,j)=1]ijd
将ddd提出
∑d=1nd∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]ij
\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]ij
d=1∑ndi=1∑⌊dn⌋j=1∑⌊dm⌋[gcd(i,j)=1]ij
莫比乌斯函数的性质∑d∣nμ(d)=[n=1]\sum_{d|n} \mu(d)=[n=1]∑d∣nμ(d)=[n=1],可以将布尔表达式转化为整除,于是即求
∑d=1nd∑i=1⌊nd⌋∑j=1⌊md⌋ij∑t∣gcd(i,j)μ(t)
\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij\sum_{t|gcd(i,j)}\mu(t)
d=1∑ndi=1∑⌊dn⌋j=1∑⌊dm⌋ijt∣gcd(i,j)∑μ(t)
再次优先枚举因子ttt,此时知道因子,只需找因子的所有倍数的和,倍数有
t,2t,3t,....
t,2t,3t,....
t,2t,3t,....
那么对于nnn,小于nnn的倍数有⌊nt⌋\lfloor\frac{n}{t}\rfloor⌊tn⌋个,求和即为
t(1+2+3+...+⌊nd⌋)
t(1+2+3+...+\lfloor\frac{n}{d}\rfloor)
t(1+2+3+...+⌊dn⌋)
令sum(i)=1+2+...+i=(1+i)(i)2sum(i)=1+2+...+i=\frac{(1+i)(i)}{2}sum(i)=1+2+...+i=2(1+i)(i),那么原式可以写作
∑d=1nd∑tmin{⌊nd⌋,⌊md⌋}t2μ(t)sum(⌊nd⌋t)sum(⌊md⌋t)
\sum_{d=1}^{n}d\sum_{t}^{min\{\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor\}}t^2\mu(t) sum(\frac{\lfloor\frac{n}{d}\rfloor}{t})sum(\frac{\lfloor\frac{m}{d}\rfloor}{t})
d=1∑ndt∑min{⌊dn⌋,⌊dm⌋}t2μ(t)sum(t⌊dn⌋)sum(t⌊dm⌋)
在题目内,可以交换n,mn,mn,m,那么不妨假设n≤mn\leq mn≤m,即求
∑d=1nd∑t⌊nd⌋t2μ(t)sum(⌊ndt⌋)sum(⌊mdt⌋)
\sum_{d=1}^{n}d\sum_{t}^{\lfloor\frac{n}{d}\rfloor}t^2\mu(t) sum(\lfloor\frac{n}{dt}\rfloor)sum(\lfloor\frac{m}{dt}\rfloor)
d=1∑ndt∑⌊dn⌋t2μ(t)sum(⌊dtn⌋)sum(⌊dtm⌋)
枚举T=dtT=dtT=dt,由于此时ttt是TTT的一个因子,此时枚举ddd为TTT的因子,地位相当于上式的ttt
∑T=1nsum(⌊nT⌋)sum(⌊mT⌋)T∑d∣Tdμ(d)
\sum_{T=1}^{n}sum(\lfloor\frac{n}{T}\rfloor)sum(\lfloor\frac{m}{T}\rfloor)T\sum_{d|T}d\mu(d)
T=1∑nsum(⌊Tn⌋)sum(⌊Tm⌋)Td∣T∑dμ(d)
其中sum(⌊nT⌋)sum(⌊mT⌋)sum(\lfloor\frac{n}{T}\rfloor)sum(\lfloor\frac{m}{T}\rfloor)sum(⌊Tn⌋)sum(⌊Tm⌋)可以整除分块解决,此时的整除分块应该保证⌊nT⌋\lfloor\frac{n}{T}\rfloor⌊Tn⌋与⌊mT⌋\lfloor\frac{m}{T}\rfloor⌊Tm⌋都在这个范围是不变的,所以rrr端点需要两次结果取小,剩下只需要快速求得
f(T)=∑d∣Tdμ(d)
f(T)=\sum_{d|T}d\mu(d)
f(T)=d∣T∑dμ(d)
一般数论函数都是积性的,所以我们考虑给一个求出来的f(T)f(T)f(T)的TTT添加一个质因子ppp
-
如果TTT含有ppp,由于当一个数nnn出现了质因子的平方的时候,μ(n)=0\mu(n)=0μ(n)=0,所以这个p,Tp,Tp,T的分解方式不能给μ(Tp)\mu(Tp)μ(Tp)带来贡献,所以f(Tp)=f(T)f(Tp)=f(T)f(Tp)=f(T)
-
如果TTT不含有ppp,由于积性函数,f(Tp)=f(T)f(p)f(Tp)=f(T)f(p)f(Tp)=f(T)f(p)
这样就可以线性筛出所有的f(T)f(T)f(T),只需要对Tf(T)Tf(T)Tf(T)做前缀和即可
设一个质因子为ppp,由于μ(p)=1−1=0\mu(p)=1-1=0μ(p)=1−1=0,于是设定f(p)=μ(p)f(p)=\mu(p)f(p)=μ(p)时即筛μ(i)\mu(i)μ(i),上述方法只是一个筛μ(i)\mu(i)μ(i)的改动版
总时间复杂度O(n)O(n)O(n)
// #include<bits/stdc++.h>
#include<iostream>
#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#include<vector>
#include<bitset>
#include<map>
using namespace std;
using ll=long long;
const int N=1e7+5,inf=0x3fffffff;
const long long INF=0x3fffffffffffffff,mod=20101009;
bool nt[N];
int prime[N],cnt;
ll n,m,ans,f[N],sum[N];
void make_prime()
{
f[1]=1;
for(int i=2;i<=n;i++)
{
if(!nt[i]) prime[++cnt]=i,f[i]=((1-i)%mod+mod)%mod;
for(int j=1;j<=cnt&&i*prime[j]<=n;j++)
{
nt[i*prime[j]]=true;
if(i%prime[j]==0)
{
f[i*prime[j]]=f[i];
break;
}
else f[i*prime[j]]=f[i]*f[prime[j]]%mod;
}
}
for(int i=1;i<=n;i++) f[i]=(f[i-1]+1ll*i*f[i]%mod)%mod;
for(int i=1;i<=m;i++) sum[i]=(sum[i-1]+i)%mod;
}
ll query(ll l,ll r){
return ((f[r]-f[l-1])%mod+mod)%mod;
}
int main()
{
cin>>n>>m;
if(n>m) swap(n,m);
make_prime();
for(ll l=1,r,val;l<=n;l=r+1)
{
val=n/l; r=val?min(n/val,n):n;
val=m/l; r=min(r,val?min(m/val,m):m);
ans=(ans+sum[n/l]*sum[m/l]%mod*query(l,r)%mod)%mod;
}
cout<<ans;
return 0;
}