BZOJ3994: [SDOI2015]约数个数和

Description

 设d(x)为x的约数个数，给定N、M，求  

 

 

Input

输入文件包含多组测试数据。

第一行，一个整数T，表示测试数据的组数。

接下来的T行，每行两个整数N、M。

 

Output

 T行，每行一个整数，表示你所求的答案。

 

Sample Input

2
7 4
5 6

Sample Output

110
121

HINT

 1<=N, M<=50000

1<=T<=50000

一开始写了个这东西，结果搞不动心态崩了，下面：

$$sum_{i=1}^{n}sum_{j=1}^{m}ij-phi(ij)$$

不知道可不可以这样做，有会的神犇请动动手指D飞我这个蒟蒻，QQ：1278795351

换做法

考虑替换d(ij)，这里有个公式：

$$d(nm)=sum_{x|n}^{n}sum_{y|m}^{m}1[gcd(x,y)==1]$$

详细证明请看钊爷博客%%%：传送门

推一波柿子：
$$sum_{i=1}^{n}sum_{j=1}^{m}sum_{x|n}^{n}sum_{y|m}^{m}1[gcd(x,y)==1]$$

$$==    sum_{x=1}^{n}sum_{y=1}^{m}frac{n}{x}*frac{m}{y}[gcd(x,y)==1]$$

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|x}mu(d)$$
提d出来：$$== sum_{d=1}^{n}mu(d) sum_{d|x}sum_{d|y}frac{n}{x}*frac{m}{y}$$

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

然后预处理后面两项，设$$S(x)=sum_{i=1}^{x}frac{x}{i}$$

 yy一下，这就是求x里面有几个i，那其实就是：

$$S(x)=sum_{i=1}^{x}d(i)$$

然后分块加速

终于写完了。。

代码如下：

//MT_LI
#include<cmath>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
using namespace std;
typedef long long ll;
const int maxn=50000;
int mu[51000],v[51000],prime[51000];
ll d[51000];
void get_Mobius()
{
    memset(v,false,sizeof(v));
    prime[0]=0;mu[1]=1;
    for(int i=2;i<=maxn;i++)
    {
        if(!v[i])
        {
            mu[i]=-1;
            prime[++prime[0]]=i;
        }
        for(int j=1;j<=prime[0]&&i*prime[j]<=maxn;j++)
        {
            v[i*prime[j]]=true;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<=maxn;i++)mu[i]+=mu[i-1];
}
int T;
int n,m;
ll get(int x)
{
    ll ans=0;
    for(int i=1;i*i<=x;i++)
    {
        if(x%i==0)
        {
            ans++;
            if(i*i!=x)ans++;
        }
    }
    return ans;
}
int main()
{
    get_Mobius();d[1]=1;
    for(int i=2;i<=maxn;i++)d[i]=(ll)get(i)+d[i-1];
    scanf("%d",&T);
    while(T--)
    {
        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)(mu[j]-mu[i-1])*d[n/i]*d[m/i];
        }
        printf("%lld
",ans);
    }
    return 0;
}