AcWing 215. 破译密码

题目传送门

参考题解：https://www.acwing.com/solution/content/49435/

一、题目分析

达达正在破解一段密码，他需要回答很多类似的问题：
对于给定的整数 \(a,b\) 和 \(d\)，有多少正整数对 \(x,y\)，满足 \(x≤a，y≤b\)，并且 \(gcd(x,y)=d\)。作为达达的同学，达达希望得到你的帮助。

我们分析一下，题目很明了，要求的就是\(gcd(x,y)=d\)的个数，转化一下：

\[\large x′=x/d,y′=y/d \]

\[\large [gcd(x,y)=d] \Leftrightarrow [gcd(x′,y′)=1] \]

注：上面的"[]"，表示的是符合这样条件的\((x,y)\)的数对个数,这两个式子的数对\((x,y)\)和\((x',y')\)是一一对应的。

这两个是等价的，但是看的时候我还是愣了一下(还是比较笨的)，但是仔细思考下，感觉还是很容易证明的。

证明：
首先，\(x , y\) 的最大公约数是\(d\)，那么我们除去最大公约数，两个数就互质了，这个毋庸置疑，对于每个\(gcd(x,y)=d\)都如此操作，那么可以证明，这两个的数量是一一对应的。
(这么显然，我还要思考一会，真滴笨！\(OVO\))

那么我们的问题转化为：\(x≤a/d\) , \(y≤b/d\) ，\(gcd(x,y)=1\)的个数，我们考虑：

补集思想:互质的数 = 总数 - 不是互质的数

总数很容易，\(a′=a/d,b′=b/d\) ， 那么总数就是\(a′?b′\)，下面分析下我们要减去不互质的数：

\[\large \sum_{i=2}^{min(a',b')}?\frac{a'}{i}?*?\frac{b'}{i}?*mobius[i] \]

式子上午看的，人是上午没的，我直接喵喵喵？？ 愣了一会，发现显然(我是笨比).

证明：
我们要减去的是什么呢？是不互质的数，不互质数是什么呢？是最大公因数不为\(1\)的数，不为\(1\)，那可以为几呢？不是\(1\)就行（废话），我们枚举最大公因数不为\(1\)的个数即可。枚举每一个可能的最大公因数，并且计算个数。
当最大公因数是只有一个质因子的时候，我们可以列出：

\[\large ?\frac{a'}{2}?* ?\frac{b'}{2}?+?\frac{a'}{3}?*?\frac{b'}{3}?+..... \]

最大公因数是\(2\)的两对数，实际上就等\(2\)的倍数个数相乘,这两个个数是等同的，同理\(3\)也是，但是我们要考虑重复的情况，如\(6\)，我们\(2，3\)，都枚举了一次，那么会重复，根据容斥原理，我们要减去，最后我们发现，其实符号就是\(Mobius\)函数（这个证明是显然的）。以此类推，枚举到\(a′\)和\(b′\)中较小的就可以了，因为大于较小的，则式子为\(0\)，我们就能推出上面的式子了。(被自己蠢哭了)

那么答案呼之欲出：

\[\large a'b'+\sum_{i=2}^{min(a',b')}?\frac{a'}{i}?*?\frac{b'}{i}?*mobius[i] \]

我们把这玩意合体一下

\[\large \sum_{i=1}^{min(a',b')}?\frac{a'}{i}?*?\frac{b'}{i}?*mobius[i] \]

（当当当当！）

然后我们发现，高兴的还是有点早，每次都是 \(O(n)\)，总体是\(O(n^2)\)，时间复杂度很大，我们忍受不了。

然后神奇的东西就出现了，我们称之为 “优雅的暴力——分块”。

下面来分析下分块，这个还好！没怎么卡。

首先，我们发现式子里面有很多取整，那么实际上数据应该是一段一段的，我们就可以考虑通过这点来加速，既然这一段都是一样的，我们直接一段一段的算，发现效率可观，可以达到\(O(n\sqrt{n})\),下面写下分块的原理，我们定义一个 \(get(x)\) 函数，这个函数和\(lowbit\)函数一样，虽然短小，但是精悍，他主要做的事是：传回当前一样的一段中，下标最大的那个,也就是\(\displaystyle \frac{n}{get(x)+1}\)就会变小。

int g(int n , int x) {
    return n / (n / x);
}

\(x\)为左边界时，右边界的计算办法：

\[\LARGE g(x)=?\frac{a}{?\frac{a}{x}?}? \]

注：数论分块的证明及相关内容

那么我们就可以愉快的写代码了～

#include 
using namespace std;

typedef long long LL;
const int N = 50010;

//线性筛法求莫比乌斯函数(枚举约数)
int mu[N], sum[N]; // 莫比乌斯函数的前缀和
int primes[N], cnt;
bool st[N];
void get_mobius(LL n) {
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!st[i]) {
            primes[cnt++] = i;
            mu[i] = -1; //奇数个质因子，现在只有1个质因子，所以函数值是-1
        }
        for (int j = 0; primes[j] * i <= n; j++) {
            int t = primes[j] * i;
            st[t] = true; //把t筛掉

            // t里有primes[j],而i里如果还有一个primes[j],那么最少有2个及以上的primes[j],根据mobius函数定义，此时函数值为0
            if (i % primes[j] == 0) {
                mu[t] = 0;
                break;
            }
            mu[t] = mu[i] * -1;
            //因为执行到这里primes[j]这个质因子只有1个，所以整个莫比乌斯函数里有没有某个质因子的个数大于1个，取决于i的质因子个数
        }
    }
    // 维护莫比乌斯函数前缀和
    for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mu[i];
}

int main() {
    //筛法求莫比乌斯函数
    get_mobius(N - 1);

    int T;
    cin >> T;

    while (T--) {
        int a, b, d;
        cin >> a >> b >> d;
        //套路啊，满满的套路，直接先用最大公约数a/gcd(a,b)=a',b/gcd(a,b)=b',映射到a',b'
        a /= d, b /= d;

        // n为 min(a', b')
        int n = min(a, b);

        LL res = 0;

        // l r, 是每一段的左右边界
        // 每次只能取较小的那个上界作为这一段的右端点r
        // 然后下次迭代时下一段的左端点就是r + 1
        for (int l = 1, r; l <= n; l = r + 1) { //分块大法
            r = min(n, min(a / (a / l), b / (b / l)));
            res += (sum[r] - sum[l - 1]) * (LL)(a / l) * (b / l);
        }
        printf("%lld\n", res);
    }
    return 0;
}