BZOJ2154: Crash的数字表格 BZOJ2693: jzptab

【传送门：BZOJ2154&BZOJ2693】

---

简要题意：

　　给出n,m，求$\sum_{i=1}^{n}\sum_{j=1}^{m}LCM(i,j)$

---

题解：

　　莫比乌斯反演（因为BZOJ2693是多组数据，数据强一点，所以代码用BZOJ2693的）

　　设n<m，原式等于$\sum_{i=1}^{n}\sum_{j=1}^{m}i*j/gcd(i,j)$

　　然后枚举d值作为i和j的gcd，得到$$\sum_{d=1}^{n}\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{i*j}{d}[gcd(i,j)==d]$$

　　因为gcd(i,j)==d，所以gcd(i/d,j/d)==1，得到$$\sum_{d=1}^{n}d*\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{m}{d}}i*j[gcd(i,j)==1]$$

　　因为莫反的性质：$\sum_{d|x}\mu(d)=[x==1]$，所以转化为$$\sum_{d=1}^{n}d*\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{m}{d}}i*j*\sum_{t|gcd(i,j)}\mu(t)$$

　　交换和式得到$$\sum_{d=1}^{n} d* \sum_{t=1}^{\frac{n}{d}} \mu(t) * \sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{m}{d}} i*j[gcd(i,j)==t]$$

　　$$\sum_{d=1}^{n}d*\sum_{t=1}^{\frac{n}{d}}t^{2}*\mu(t)*\sum_{i=1}^{\frac{n}{dt}}\sum_{j=1}^{\frac{m}{dt}}i*j$$

　　设$T=dt$，$S(x)=\sum_{i=1}^{x}i$，得到$$\sum_{d=1}^{n}\sum_{t=1}^{\frac{n}{d}}T*t*\mu(t)*S(\frac{n}{T})*S(\frac{m}{T})$$

　　将$S(\frac{n}{T})*S(\frac{m}{T})$提前，得到$$\sum_{T=1}^{n}T*S(\frac{n}{T})*S(\frac{m}{T})\sum_{d|T}d*\mu(d)$$

　　因为$S(\frac{n}{T})*S(\frac{m}{T})$可以前缀和预处理，显然我们只要将$\sum_{d|T}d*\mu(d)$快速求出就可以了

　　设$F(T)=\sum_{d|T}d*\mu(d)$，显然是一个积性函数，在线性筛的时候求就行了

　　然后将$T*F(T)$求前缀和，然后整除分块加速就能过了

---

参考代码：

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<cmath>
#include<algorithm>
using namespace std;
typedef long long LL;
int prime[11000000],v[11000000];
LL f[11000000],sum[11000000];
LL Mod=1e8+9;
void pre(int n)
{f[1]=1;int tot=0;for(int i=2;i<=n;i++){if(v[i]==0){v[i]=i;prime[++tot]=i;f[i]=(1-i+Mod)%Mod;}for(int j=1;j<=tot;j++){if(prime[j]>v[i]||prime[j]>n/i) break;v[i*prime[j]]=prime[j];if(i%prime[j]==0){f[i*prime[j]]=f[i]%Mod;break;}else f[i*prime[j]]=f[i]*f[prime[j]]%Mod;}}for(int i=1;i<=n;i++) f[i]=(f[i]*LL(i)%Mod+f[i-1])%Mod;for(int i=1;i<=n;i++) sum[i]=(sum[i-1]+LL(i))%Mod;
}
int main()
{pre(10000000);int T;scanf("%d",&T);while(T--){int n,m;scanf("%d%d",&n,&m);if(n>m) swap(n,m);LL ans=0;for(int i=1,j;i<=n;i=j+1){j=min(n/(n/i),m/(m/i));ans=(ans+(f[j]-f[i-1]+Mod)%Mod*sum[n/i]%Mod*sum[m/i]%Mod)%Mod;}printf("%lld\n",ans);}return 0;
}