写一篇更适合蒟蒻体质的题解

题目简述：求区间 $[a,b]$ 中与 $n$ 互质的数的个数

如果要求的数是区间 $[1,n]$ 中与 $n$ 互质的数，用欧拉函数就能很简单地解决，但是区间的范围显然不是 $[1,n]$，所以我们要求的其实是 $[1,m]$ 中与 $n$ 互质的数的个数这样一个问题。

求解这个问题的话我们可以对 $n$ 进行质因数分解

• 当 $n = p^k$ 的形式的时候，$[1,m]$ 中与 $n$ 互质的数字也一定不是 $p$ 的倍数，所以 $[1,m]$ 中与 $n$ 互质的数就是 $m$ 减去 $p$ 的倍数的个数，也就是 $m - {m \over p}$
• 再看到 $n = {p_1}^{k_1}· {p_2}^{k_2}$ 的形式的时候，此时可以发现 ${m \over p}$ 就变成了 $p_1$ 的倍数个数加上 $p_2$ 的倍数个数再减去他们被多加了一次的 $p_1$ 和 $p_2$ 的公倍数，也就是${m \over {p_1}}+{m \over {p_2}}-{m \over {p_2p_1}}$ ，可以发现这是一个容斥原理的问题
• 所以如果把这个式子推广到 $n = {p_1}^{k_1}· {p_2}^{k_2}...{p_x}^{k_x}$ 的时候， ${m \over p}$ 就是 ${m \over {p_i}}+{m \over {p_j}}+{m \over {p_k}}-{m \over {p_ip_j}}-{m \over {p_jp_k}}-{m \over {p_ip_k}}+{m \over {p_ip_jp_k}}$ 的形式，通过枚举这个集合并对他进行容斥原理奇数加偶减的操作，就可以求出 ${m \over p}$，这个结果就是 $[1,m]$ 中不与 $n$ 互质的数的个数，最后用 $m$ 减去它，就是 $[1,m]$ 中与 $n$ 互质的数的个数。

对于区间 $[a,b]$ 直接通过上述方式算出 $[1,b]$ 以及 $[1,a - 1]$ 中不与 $n$ 互质的数的个数 ${m \over {p_b}}，{m \over {p_a}}$，然后再用 $[a,b]$ 的元素个数 $(b-a+1)$ 减去 ${m \over {p_b}}$ 再加上多减了的 ${m \over {p_a}}$ 就是答案

具体的

1.先筛出一定数量的质数（最好是 $\le \sqrt{N}$），再对n进行分解质因数

bool vis[MAXN];
vector<int> Prime;
void get_prime(int n) {
	for (int i = 2; i <= n; i++) {
		if (!vis[i]) Prime.push_back(i);
		for (auto p : Prime) {
			if (p * i > n) break;
			vis[i * p] = true;
			if (i % p == 0) break;
		}
	}
}
LL num, factor[MAXN];
void get_factor(int n) {
	num = 0;
	for (auto p : Prime) {
		if (p * p == n) break;
		if (n % p != 0) continue;
		factor[num ++] = p;
		while (n % p == 0)
			n /= p;
	}
	if (n != 1) factor[num ++] = n;
}

2. 考虑到直接用DFS枚举子集复杂度太高，所以这里使用二进制枚举子集

LL get_mult(LL m, LL num) {
	LL res = 0;
	for (int i = 1; i < (1 << num); i++) {
		LL tot = 0;
		LL tmp = 1;
		for (int j = 0; j < num; j++) 
			if (i & (1 << j)) 
				tot ++, tmp *= factor[j];
		if (tot & 1) res += m / tmp;
		else res -= m / tmp;
	}
	return res;
}

3. 处理输出

void solve() {
	LL a, b, n;
	scanf ("%lld%lld%lld", &a, &b, &n);	
	get_factor(n);
	++ case_num;
	LL ans = (b - a + 1);
	ans -= get_mult(b, num) - get_mult(a - 1, num);
	printf ("Case #%d: %lld", case_num, ans);
}

完整代码

#include <iostream>
#include <vector>
#include <algorithm>
#define MAXN (100010)

using namespace std;
typedef long long LL;
int T, case_num;

bool vis[MAXN];
vector<int> Prime;
void get_prime(int n) {
	for (int i = 2; i <= n; i++) {
		if (!vis[i]) Prime.push_back(i);
		for (auto p : Prime) {
			if (p * i > n) break;
			vis[i * p] = true;
			if (i % p == 0) break;
		}
	}
}

LL num, factor[MAXN];
void get_factor(int n) {
	num = 0;
	for (auto p : Prime) {
		if (p * p == n) break;
		if (n % p != 0) continue;
		factor[num ++] = p;
		while (n % p == 0)
			n /= p;
	}
	if (n != 1) factor[num ++] = n;
}

LL get_mult(LL m, LL num) {
	LL res = 0;
	for (int i = 1; i < (1 << num); i++) {
		LL tot = 0;
		LL tmp = 1;
		for (int j = 0; j < num; j++) 
			if (i & (1 << j)) 
				tot ++, tmp *= factor[j];
		if (tot & 1) res += m / tmp;
		else res -= m / tmp;
	}
	return res;
}

void solve() {
	LL a, b, n;
	scanf ("%lld%lld%lld", &a, &b, &n);	
	get_factor(n);
	++ case_num;
	LL ans = (b - a + 1);
	ans -= get_mult(b, num) - get_mult(a - 1, num);
	printf ("Case #%d: %lld", case_num, ans);
}

int main() {
	get_prime(MAXN - 1);
	scanf ("%d", &T);
	while (T --) solve();
	return 0;
}

PS.POJ的和HDU的在这里居然都提交不了，太可恶，有人知道怎么提交嘛