BZOJ3561 DZY Loves Math VI

原题链接：https://www.lydsy.com/JudgeOnline/problem.php?id=3561

DZY Loves Math VI

Description

给定正整数n,m。求

∑i=1n∑j=1mlcm(i,j)gcd(i,j)∑i=1n∑j=1mlcm(i,j)gcd(i,j)

\sum_{i=1}^{n}\sum_{j=1}^{m}lcm(i,j)^{gcd(i,j)}

Input

一行两个整数n,m。

Output

一个整数，为答案模1000000007后的值。

Sample Input

5 4

Sample Output

424

HINT

数据规模：

1<=n,m<=500000，共有3组数据。

题解

∑i=1n∑j=1mlcm(i,j)gcd(i,j)=∑i=1n∑j=1m(ijgcd(i,j))gcd(i,j)=∑d=1min(n,m)∑i=1n∑j=1m(ijd)d[d=gcd(i,j)]=∑d=1min(n,m)∑i=1⌊nd⌋∑j=1⌊md⌋(d2ijd)d[1=gcd(i,j)]=∑d=1min(n,m)dd∑i=1⌊nd⌋∑j=1⌊md⌋(ij)d∑d′|gcd(i,j)μ(d′)=∑d=1min(n,m)dd∑d′=1⌊min(n,m)d⌋μ(d′)∑i=1⌊nd⌋∑j=1⌊md⌋(ij)d[d′|gcd(i,j)]=∑d=1min(n,m)dd∑d′=1⌊min(n,m)d⌋μ(d′)d′2d∑i=1⌊ndd′⌋id∑j=1⌊mdd′⌋jd∑i=1n∑j=1mlcm(i,j)gcd(i,j)=∑i=1n∑j=1m(ijgcd(i,j))gcd(i,j)=∑d=1min(n,m)∑i=1n∑j=1m(ijd)d[d=gcd(i,j)]=∑d=1min(n,m)∑i=1⌊nd⌋∑j=1⌊md⌋(d2ijd)d[1=gcd(i,j)]=∑d=1min(n,m)dd∑i=1⌊nd⌋∑j=1⌊md⌋(ij)d∑d′|gcd(i,j)μ(d′)=∑d=1min(n,m)dd∑d′=1⌊min(n,m)d⌋μ(d′)∑i=1⌊nd⌋∑j=1⌊md⌋(ij)d[d′|gcd(i,j)]=∑d=1min(n,m)dd∑d′=1⌊min(n,m)d⌋μ(d′)d′2d∑i=1⌊ndd′⌋id∑j=1⌊mdd′⌋jd

\begin{align*} &\sum_{i=1}^{n}\sum_{j=1}^{m}lcm(i,j)^{gcd(i,j)}\\ &=\sum_{i=1}^{n}\sum_{j=1}^{m}(\frac{ij}{gcd(i,j)})^{gcd(i,j)}\\ &=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{n}\sum_{j=1}^{m}(\frac{ij}{d})^{d}[d=gcd(i,j)]\\ &=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}(\frac{d^2ij}{d})^{d}[1=gcd(i,j)]\\ &=\sum_{d=1}^{min(n,m)}d^d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}(ij)^{d}\sum_{d'|gcd(i,j)}\mu(d')\\ &=\sum_{d=1}^{min(n,m)}d^d\sum_{d'=1}^{\lfloor\frac{min(n,m)}{d}\rfloor}\mu(d')\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}(ij)^{d}[d'|gcd(i,j)]\\ &=\sum_{d=1}^{min(n,m)}d^d\sum_{d'=1}^{\lfloor\frac{min(n,m)}{d}\rfloor}\mu(d')d'^{2d}\sum_{i=1}^{\lfloor\frac{n}{dd'}\rfloor}i^d\sum_{j=1}^{\lfloor\frac{m}{dd'}\rfloor}j^{d}\\ \end{align*}

因为dd′dd′dd'后面还跟了一个ididi^d之类的奇奇怪怪的东西，所以用TTT来替换并无卵用，但是数据范围只有5e4" role="presentation" style="position: relative;">5e45e45e4啊，如果直接枚举的话复杂度为n1+n2+n3+⋯+nn−1+nn≈n ln(n)n1+n2+n3+⋯+nn−1+nn≈nln(n)\frac{n}{1}+\frac{n}{2}+\frac{n}{3}+\cdots+\frac{n}{n-1}+\frac{n}{n}\approx n\ ln(n)，稳了稳了。

代码

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int M=5e5+5,mod=1000000007;
int p[M/3],miu[M],n,m,nd,md,dd;
ll ans,tmp,sum[M],val[M];
bool isp[M];
ll power(ll x,ll p){ll r=1;for(;p;p>>=1,x=x*x%mod)if(p&1)r=r*x%mod;return r;}
void pre()
{miu[1]=isp[1]=1;for(int i=2,t;i<M;++i){if(!isp[i])p[++p[0]]=i,miu[i]=-1;for(int j=1;j<=p[0],i*p[j]<M;++j){isp[i*p[j]]=1;if(i%p[j]==0)break;miu[i*p[j]]=-miu[i];}}
}
void in(){scanf("%d%d",&n,&m);}
void ac()
{if(n>m)swap(n,m);pre();for(int i=1;i<=m;++i)val[i]=1;for(int i=1,j;i<=n;++i){nd=n/i,md=m/i,tmp=0;for(j=1;j<=md;++j)val[j]=val[j]*j%mod,sum[j]=(sum[j-1]+val[j])%mod;for(j=1;j<=nd;++j)(tmp+=miu[j]*val[j]%mod*val[j]%mod*sum[nd/j]%mod*sum[md/j]%mod)%=mod;(ans+=tmp*power(i,i)%mod)%=mod;}printf("%lld",(ans+mod)%mod);
}
int main(){in();ac();}