bzoj2820&&COGS2165 YY的GCD

bzoj传送门（权限）
COGS传送门
Description
神犇YY虐完数论后给傻×kAc出了一题给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少对kAc这种傻×必然不会了，于是向你来请教……多组输入
Input
第一行一个整数T 表述数据组数接下来T行，每行两个正整数，表示N, M
Output
T行，每行一个整数表示第i组数据的结果
Sample Input
2
10 10
100 100
Sample Output
30
2791
HINT
T = 10000
N, M <= 10000000
Source

题解

mobius反演。。
正经地写一篇关于反演的题解
设F(x)为n|gcd(i,j)(i,j在给定范围内)的i,j的个数
f(x)为n=gcd(i,j)的i,j的个数
根据反演，得

$$f(n)=\sum_{n|d}μ(\frac{d}{n})F(d)$$
因为以 $p[i]$(p是素数数组)为gcd等价于同除以 $p[i]$后互质
$$=\sum_{i-1}^{tot}\sum_{p[i]|d}μ(d)\lfloor\frac{n}{p[i]d}\rfloor\lfloor\frac{m}{p[i]d}\rfloor$$
令 $p[i]d=T$则
$$ans=\sum_{i=1}^{tot}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\sum_{p[i]|T}μ(\frac{T}{p[i]})$$
//其实我更喜欢把最后的 $\frac{T}{p[i]}$看作 $d$，这样在最后写代码的时候理解起来方便一些
这样以后如果我们能预处理$\sum_{p[i]|T}μ(d)$就能完成要求
我们只需要枚举每个质数，更新其所有倍数就好了
剩下的和原先的反演几乎没有区别

CODE：

#include<cstdio>
#define N 10000005
#define maxn 700000
typedef long long ll;
int prime[maxn];
bool Prime[N];
int mu[N];
ll f[N];
int T,n,m,tot;
inline int min(const int &a,const int &b){return a<b?a:b;}
inline void swap(int &a,int &b){a^=b,b^=a,a^=b;}
inline void mobius()
{
    mu[1]=1;
    for(int i=2;i<=N-5;i++)
    {
        if(!Prime[i]) prime[++tot]=i,mu[i]=-1;
        for(int j=1;j<=tot&&i*prime[j]<=N-5;j++)
        {
            Prime[i*prime[j]]=1;
            if(i%prime[j]==0) break;
            mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<=tot;i++)
      for(int j=1;j*prime[i]<=N-5;j++)
        f[j*prime[i]]+=mu[j];
    for(int i=2;i<=N-5;i++)
      f[i]+=f[i-1];
}
inline ll ask(int n,int m)
{
    if(n>m) swap(n,m);
    ll ans=0;
    for(int i=1,last;i<=n;i=last+1)
    {
        last=min(n/(n/i),m/(m/i));
        ans+=(f[last]-f[i-1])*(n/i)*(m/i);
    }
    return ans;
}
int main()
{
    scanf("%d",&T);
    mobius();
    while(T--)
    {
        scanf("%d%d",&n,&m);
        printf("%lld\n",ask(n,m));
    }
    return 0;
}