[SDOI2015]约数个数和 莫比乌斯反演

～～～题面～～～

题解：

为什么SDOI这么喜欢莫比乌斯反演，，，

首先有一个结论$$d(ij) = \sum_{x|i}\sum_{y|j}[gcd(x, y) == 1]$$
为什么呢？
首先，可以看做从两个数中分别取一些不重叠的质数的$k_{i}$次方，组成新数的方案数。
那如果有需要重叠的部分怎么办？
可以看做全都在第一个or第二个数中取。
但是一个数的次数不够怎么办呢？
相当于以1为媒介，可以简介统计到这些情况
比如$2^{3} \cdot 2^{2}  = 2^{5}$可以看成$1, 2^{3}$ + $1, 2^{2}$,这样就可以做到分别枚举了指数从1到5的情况了；
$$ans = \sum_{i = 1}^{n}\sum_{j = 1}^{m}\sum_{x | i}\sum_{y | j}[gcd(x, y) == 1]$$
$$ans = \sum_{x = 1}^{n}\sum_{y = 1}^{m} \lfloor{\frac{n}{x}}\rfloor \lfloor{\frac{m}{y}}\rfloor [gcd(x, y) == 1]$$<---改成枚举x和y，然后统计有多少对i,j中分别含有因数x，y
设$$f(d) = \sum_{x = 1}^{n}\sum_{y = 1}^{m} \lfloor{\frac{n}{x}}\rfloor \lfloor{\frac{m}{y}}\rfloor [gcd(x, y) == d]$$, $$g(x) = \sum_{x | d}^{min(n, m)}{f(d)}$$
那么$ans = f(1)$
反演一下得到:$$f(n)  = \sum_{n | d}\mu(\frac{d}{n}) \cdot g(d)$$
接下来就是怎么求$g(d)$了。
$$g(d) = \sum_{x = 1}^{n}\sum_{y = 1}^{m} \lfloor{\frac{n}{x}}\rfloor \lfloor{\frac{m}{y}}\rfloor [d | gcd(x, y)]$$
$d | gcd(x, y)$ ---> $d | x + d | y$ ,因为$x | n + y | m$ ---> $d | n + d | m$
所以枚举$d_{i}$就符合条件了，$\lfloor{\frac{n}{dx}}\rfloor$ ---> $dx$倍数的个数，即有多少i，j的搭配是可以满足$[d|gcd(x, y)]$的
$$g(d) = \sum_{x = 1}^{ \lfloor{\frac{n}{d}}\rfloor }\sum_{y = 1}^{ \lfloor{\frac{m}{d}}\rfloor } \lfloor{\frac{n}{dx}}\rfloor \lfloor{\frac{m}{dy}}\rfloor$$
于是现在就是要求$f(1)$
$$f(1) = \sum_{d = 1}^{min(n, m)}{\mu(\frac{d}{1})g(d)}$$
$$= \sum_{d = 1}^{min(n, m)}{\mu(d)g(d)}$$
$$= \sum_{d = 1}^{min(n, m)}{\mu(d) \cdot \sum_{x = 1}^{ \lfloor{\frac{n}{d}}\rfloor }\sum_{y = 1}^{ \lfloor{\frac{m}{d}}\rfloor } \lfloor{\frac{n}{dx}}\rfloor \lfloor{\frac{m}{dy}}\rfloor}$$
这个式子可以整数分块。
设$N = \frac{n}{d}$， $M = \frac{m}{d}$,那么
$$= \sum_{d = 1}^{min(n, m)}{\mu(d) \cdot \sum_{x = 1}^{ N }\sum_{y = 1}^{ M } \lfloor{\frac{N}{x}}\rfloor \lfloor{\frac{M}{y}}\rfloor}$$
$$= \sum_{d = 1}^{min(n, m)}{\mu(d) \cdot \sum_{x = 1}^{ N}{\lfloor{\frac{N}{x}}\rfloor} \sum_{y = 1}^{M}  {\lfloor{\frac{M}{y}}\rfloor}}$$
而$N = \lfloor{\frac{n}{d}}\rfloor$， $M = \lfloor{\frac{m}{d}}\rfloor$有很多相同的区间段,所以可以整数分块求。

#include<bits/stdc++.h>
using namespace std;
#define R register int
#define AC 50100
#define ac 50000
#define LL long long
int n, m, tot, t;
LL ans;
int prime[AC];
LL mu[AC], s[AC];
bool z[AC];

inline int read()
{
    int x = 0;char c = getchar();
    while(c > '9' || c < '0') c = getchar();//error!!!又一次打错了读入优化。。。是||啊
    while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
    return x;
}

void pre()
{
    int x;
    mu[1] = 1;
    for(R i = 2; i <= ac; i++)
    {
        if(!z[i]) prime[++tot] = i, mu[i] = -1;
        for(R j = 1; j <= tot; j++)
        {
            x = prime[j];
            if(x * i > ac) break;
            z[x * i] = true;
            if(!(i % x)) break;
            mu[i * x] = -mu[i]; 
        }
    }
    for(R i = 1; i <= ac; i++) mu[i] += mu[i-1];    
    int pos;
    for(R i = 1; i <= ac; i++) 
    {
        for(R j = 1; j <= i; j = pos + 1)
        {
            pos = i / (i / j);
            s[i] += (LL)(i / j) * (LL)(pos - j + 1);
        }
    }
//    for(R i = 1; i <= 30; i++) printf("%lld\n", s[i]);
}

void work()
{
    t = read();
    while(t--)
    {
        n = read(), m = read(), ans = 0;
        //scanf("%d%d", &n, &m);ans = 0;
        int pos, b = min(n, m);
        for(R i = 1; i <= b; i = pos + 1)
        {
            pos = min(n / (n / i), m / (m / i));
            ans += (mu[pos] - mu[i - 1]) * s[n / i] * s[m / i];
        }
        printf("%lld\n", ans);
    }
}

int main()
{
//    freopen("in.in", "r", stdin);
    pre();
    work();
//    fclose(stdin);
    return 0;
}