推式子太快乐啦!虽然我好蠢而且dummy和maomao好巨(划掉)
思路
莫比乌斯反演的题目
首先这题有\(O(\sqrt n)\)的做法但是我没写咕咕咕
然后就是爆推一波式子
\[
\sum_{i=1}^{n}\sum_{j=1}^{m}lcm(i,j)
\]
\[
\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{i\times j}{gcd(i,j)}
\]
设$ gcd(i,j)=d$
\[
\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\times[gcd(i,j)=1]
\]
\[
p_1=\lfloor\frac{n}{d}\rfloor\\ p_2=\lfloor\frac{m}{d}\rfloor
\]
设\(f(x)\)为满足条件\(1 \le i \le p_1\)且$1 \le j \le p_2 \(且\)[gcd(i,j)=x]\(的\)i\times j$的和
则可以推出\(F(x)\)为
\[
F(x)=\sum_{x|d}f(d)
\]
所以\(F(x)\)为满足条件\(1 \le i \le p_1\)且$1 \le j \le p_2 \(且\)[x|gcd(i,j)]\(的\)i\times j$的和
所以
\[
F(x)=x^2\times \sum_{i=1}^{\lfloor\frac{p_1}{x}\rfloor}i\times \sum_{j=1}^{\lfloor\frac{p_2}{x}\rfloor}j
\]
因为
\[
f(x)=\sum_{x|d}\mu(\frac{d}{x})F(d)
\]
所以
\[
f(x)=\sum_{x|d}\mu(\frac{d}{x})\times d^2 \times \sum_{i=1}^{\lfloor\frac{p_1}{d}\rfloor}i\times \sum_{j=1}^{\lfloor\frac{p_2}{d}\rfloor}j
\]
所以
\[
f(x)=\sum_{t=1}\mu(t)\times (tx)^2 \times \sum_{i=1}^{\lfloor\frac{p_1}{tx}\rfloor}i\times \sum_{j=1}^{\lfloor\frac{p_2}{tx}\rfloor}j
\]
\[
ans=\sum_{d=1}^{n}d\sum_{t=1}^{\lfloor\frac{n}{d}\rfloor}\mu(t)\times (t)^2 \times \sum_{i=1}^{\lfloor\frac{n}{td}\rfloor}i\times \sum_{j=1}^{\lfloor\frac{m}{td}\rfloor}j
\]
然后两个整除分块搞定,复杂度\(O(n)\)
然后dummy教了我一种新的计算方式
\[
\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\times[gcd(i,j)=1]
\]
发现[gcd(i,j)=1]的形式有点眼熟
直接把\(\mu\)带进去算
\[
\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\times \sum_{p|gcd(i,j)}\mu(p)
\]
算得
\[
\sum_{d=1}^nd\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}p^2\mu(p)\sum_{i=1}^{\lfloor\frac{n}{dp}\rfloor}i\sum_{j=1}^{\lfloor\frac{n}{dp}\rfloor}j
\]
代码
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;
const int INV = 10050505;
const int MOD = 20101009;
int isprime[10010000],iprime[10010000],cnt,mu[10010000],summu[10010000],n,m;
void prime(int n){
isprime[n]=true;
mu[1]=1;
for(int i=2;i<=n;i++){
if(!isprime[i])
iprime[++cnt]=i,mu[i]=-1;
for(int j=1;j<=cnt&&iprime[j]*i<=n;j++){
isprime[iprime[j]*i]=true;
mu[iprime[j]*i]=-mu[i];
if(i%iprime[j]==0){
mu[iprime[j]*i]=0;
break;
}
}
}
for(int i=1;i<=n;i++)
summu[i]=(summu[i-1]%MOD+1LL*(mu[i]%MOD+MOD)%MOD*i%MOD*i%MOD)%MOD;
}
long long calc2(int n,int m){
long long ans=0;
for(int l=1,r;l<=min(n,m);l=r+1){
r=min(n/(n/l),m/(m/l));
ans=(ans%MOD+1LL*(summu[r]-summu[l-1]%MOD+MOD)%MOD*(1+n/l)%MOD*(n/l)%MOD*INV%MOD*(1+m/l)%MOD*(m/l)%MOD*INV%MOD)%MOD;
}
return ans;
}
long long calc1(int n,int m){
long long ans=0;
for(int l=1,r;l<=min(n,m);l=r+1){
r=min(n/(n/l),m/(m/l));
ans=(ans%MOD+1LL*(r-l+1)%MOD*(r+l)%MOD*INV%MOD*calc2(n/l,m/l)%MOD)%MOD;
}
return ans;
}
int main(){
prime(10001000);
scanf("%d %d",&n,&m);
if(n<m)
swap(n,m);
long long ans=calc1(n,m);
printf("%lld",ans);
}
原文地址:https://www.cnblogs.com/dreagonm/p/10073395.html