题目地址：HDUOJ：Problem – 1695

题目描述

Given 5 integers: a, b, c, d, k, you’re to find x in a…b, y in c…d that GCD(x, y) = k. GCD(x, y) means the greatest common divisor of x and y. Since the number of choices may be very large, you’re only required to output the total number of different number pairs.
Please notice that, (x=5, y=7) and (x=7, y=5) are considered to be the same.
Yoiu can assume that a = c = 1 in all test cases.
求[1, b]中的整数x与[1, d]中的整数y使得gcd(x, y)=k的数对(x, y)的数量。注意(5, 7)和(7, 5)只能算一对。

输入输出格式

输入格式：
The input consists of several test cases. The first line of the input is the number of the cases. There are no more than 3,000 cases.
Each case contains five integers: a, b, c, d, k, 0 < a <= b <= 100,000, 0 < c <= d <= 100,000, 0 <= k <= 100,000, as described above.

输出格式：
For each test case, print the number of choices. Use the format in the example.

输入输出样例

输入样例#1：

2
1 3 1 5 1
1 11014 1 14409 9

输出样例#1：

Case 1: 9
Case 2: 736427

题解

本题需要的数学姿势：莫比乌斯反演原理及应用 | KSkun’s Blog
由于gcd(xk, yk)=k等价于gcd(x, y)=1，我们只需要求出[1, b/k]和[1, d/k]中互质的数对数即可。这样可以先行缩小b和d的范围。
我们构造函数f(n)为gcd(x, y)=n的数对数，而答案就是f(1)。这个函数求值很不好办，我们考虑构造另外一个函数F(n)为n|gcd(x, y)的数对数。我们先令p=b/k，q=d/k，很容易知道
F(n) = \frac{p}{n} \cdot \frac{q}{n}
然后研究f(n)和F(n)的关系，发现满足下面这个关系
F(n) = \sum_{n|d} f(d)
然后就可以用莫比乌斯反演
f(n) = \sum_{n|d} \mu(\frac{d}{n}) F(d)
由于我们求的是f(1)，所以实际上是这样求的
f(1) = \sum_{i=1}^{\min(p, q)} \mu(i) \cdot \frac{p}{i} \cdot \frac{q}{i}
而由于题目要求去重，我们发现只有当x和y都处于 [1, \min(b, d)] 中时才会发生重复，因此我们令p=q=\min(b/k, d/k)再求一遍f(1)，减去它的一半即可。

代码

// Code by KSkun, 2018/4
#include <cstdio>

#include <algorithm>

typedef long long LL;

inline char fgc() {
    static char buf[100000], *p1 = buf, *p2 = buf;
    return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF : *p1++;
}

inline LL readint() {
    register LL res = 0, neg = 1;
    register char c = fgc();
    while(c < '0' || c > '9') {
        if(c == '-') neg = -1;
        c = fgc();
    }
    while(c >= '0' && c <= '9') {
        res = res * 10 + c - '0';
        c = fgc();
    }
    return res * neg;
}

const int MAXN = 100005;

bool notprime[MAXN];
int prime[MAXN], miu[MAXN], ptot;

inline void sieve() {
    miu[1] = 1;
    for(int i = 2; i <= MAXN; i++) {
        if(!notprime[i]) {
            prime[ptot++] = i;
            miu[i] = -1;
        }
        for(int j = 0; j < ptot && prime[j] * i <= MAXN; j++) {
            int k = prime[j] * i;
            notprime[k] = true;
            if(i % prime[j] == 0) {
                miu[k] = 0; break;
            } else {
                miu[k] = -miu[i];
            }
        }
    }
}

int T, a, b, c, d, k;

int main() {
    sieve();
    T = readint();
    for(int kase = 1; kase <= T; kase++) {
        LL all = 0, rep = 0;
        a = readint(); b = readint(); c = readint(); d = readint(); k = readint();
        if(!k) {
            printf("Case %d: 0\n", kase);
            continue;
        }
        b /= k; d /= k;
        if(b > d) std::swap(b, d);
        for(int i = 1; i <= b; i++) {
            all += 1ll * miu[i] * (b / i) * (d / i);
        }
        for(int i = 1; i <= b; i++) {
            rep += 1ll * miu[i] * (b / i) * (b / i);
        }
        printf("Case %d: %lld\n", kase, all - rep / 2);
    }
    return 0;
}