BZOJ3994: [SDOI2015]约数个数和

【传送门：BZOJ3994】

---

简要题意：

　　给出n,m，设d(x)为x的约数个数，求$sum_{i=1}^{n}sum_{j=1}^{m}d(i*j)$

---

题解：

　　莫比乌斯反演，设n<m

　　yy一下可以发现$d(i*j)=sum_{x|i}sum_{y|j}1[gcd(x,y)==1]$

　　然后原式就等于$$sum_{i=1}^{n}sum_{j=1}^{m}sum_{x|i}sum_{y|j}1[gcd(x,y)==1]$$

　　先枚举x,y，得到$$sum_{x=1}^{n}sum_{y=1}^{m}frac{n}{x}*frac{m}{y}[gcd(x,y)==1]$$

　　因为$sum_{d|x}mu(d)=[x==1]$

　　所以得到$$sum_{x=1}^{n}sum_{y=1}^{m}frac{n}{x}*frac{m}{y}sum_{d|gcd(x,y)}mu(d)$$

　　先枚举d，得到$$sum_{d=1}^{n}mu(d)sum_{x=1}^{n}[d|x]sum_{y=1}^{m}[d|y]frac{n}{x}*frac{m}{y}$$

　　$$sum_{d=1}^{n}mu(d)sum_{x=1}^{frac{n}{d}}frac{n}{d*x}*sum_{y=1}^{frac{m}{d}}frac{m}{d*y}$$

　　显然最前面一项可以预处理，设$S(x)=sum_{i=1}^{x}frac{x}{i}$，得到$$sum_{d=1}^{n}mu(d)*S(frac{n}{d})*S(frac{m}{d})$$

　　S(x)可以$O(n*sqrt{n})$预处理，但是实际上可以在线性筛的时候顺便处理

　　因为我们发现$S(x)=sum_{i=1}^{x}d(i)$，只要求出d(i)就可以了

---

参考代码：

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<cmath>
#include<algorithm>
using namespace std;
typedef long long LL;
int miu[51000],prime[51000],v[51000];
LL d[51000];int c[51000];
void pre(int n)
{
    int m=0;miu[1]=d[1]=c[1]=1;
    for(int i=2;i<=n;i++)
    {
        if(v[i]==0)
        {
            v[i]=i;c[i]=1;d[i]=2;
            prime[++m]=i;
            miu[i]=-1;
        }
        for(int j=1;j<=m;j++)
        {
            if(prime[j]>v[i]||prime[j]>n/i) continue;
            v[i*prime[j]]=prime[j];
            if(i%prime[j]==0) miu[i*prime[j]]=0,d[i*prime[j]]=d[i]/(c[i]+1)*(c[i]+2),c[i*prime[j]]=c[i]+1;
            else miu[i*prime[j]]=-miu[i],d[i*prime[j]]=d[i]*d[prime[j]],c[i*prime[j]]=1;
        }
    }
    for(int i=1;i<=n;i++) miu[i]+=miu[i-1],d[i]+=d[i-1];
}
int main()
{
    pre(50000);
    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+=(LL)(miu[j]-miu[i-1])*d[n/i]*d[m/i];
        }
        printf("%lld
",ans);
    }
    return 0;
}