【多解 · 代码超详解】POJ 3292 Semi-prime H-numbers

本文介绍了一种特殊数列——H数列，以及在此数列中定义的H质数与H半质数概念。通过两种算法详细阐述了如何在H数列中寻找H半质数，并提供了具体代码实现，包括打质数表模板的运用及H半质数的高效计数。

一、题目描述

数列 {4n+1}，n = 0，1，2，……叫做 H 数，它们是：1，5，9，13，17，21，……。显然 H 数列对乘运算封闭。

集合对运算封闭：见《离散数学 第2版》（屈婉玲、耿素云、张立昂）9.1 二元运算及其性质。
知识梳理：https://blog.csdn.net/COFACTOR/article/details/100145709

在 H 数的范围内定义 H 质数：只有 1 和其本身能将自己整除的数。
1 是唯一的单位元。剩下的数是 H 合数。
H 半质数是两个 H 质数（相同或不同）的积。
输入若干行，每行一个 H 数 h（h ≤ 1 000 001），问 1 到 h 共有多少个 H 半质数。
输入 0 ，程序直接结束。
样例：
输入

输出

21 0
85 5
789 62

二、算法分析说明与代码编写指导

将打质数表的模板直接套用即可，不过范围从每个正整数变成了数列 {4n + 1}。
打好表后：
法一对每个 H 数 h 都尝试分解，分解到 h = 1 或得到了两个 H 质数时即刻结束。
法二枚举两个 H 质数的积，并将 bitset 中对应的位置标记 1，代表该积为 H 半质数。当积大于 1 000 001 时，跳过本轮循环。最后遍历一次 bitset 统计 1 到每个数的 H 半质数的数量。

三、AC 代码

法一：（172 ms）

#include<cstdio>
#include<algorithm>
#include<cmath>
#include<bitset>
#pragma warning(disable:4996)
using namespace std;
unsigned HPrime[20000] = { 5,9 }, _PTy, MaxHPrime, * HPrime_end = HPrime + sizeof(HPrime) / sizeof(HPrime[0]);
unsigned x, num[1000002], n = 0; bitset<250000> isHPrime;
inline void genprime() {
	unsigned a = 13, t; bool flag = true;
	for (unsigned* i = HPrime + 2; i != HPrime_end;) {
		t = sqrt(a); flag = true;
		for (unsigned* j = HPrime; *j <= t; ++j)
			if (a % *j == 0) { flag = false; break; }
		if (flag) { *i = a, ++i, isHPrime[a] = 1; }
		a += 4;
	}
	MaxHPrime = *(HPrime_end - 1);
}
template<typename _Ty> inline bool pfact(_Ty x) {
	static _Ty t, L = MaxHPrime * MaxHPrime;
	t = min((_Ty)sqrt(x), MaxHPrime);
	for (unsigned* d = HPrime; *d <= t; ++d) {
		if (x % *d == 0) {
			x /= *d; return isHPrime[x];
		}
	}
	return false;
}
int main() {
	genprime(); isHPrime[5] = isHPrime[9] = 1;
	for (unsigned i = 5; i <= 1000001; i += 4) {
		if (pfact(i) == true) { ++n; }
		num[i] = n;
	}
	for (;;) {
		scanf("%u", &x); if (x == 0)return 0;
		printf("%u %u\n", x, num[x]);
	}
}

法二：（16 ms）

#include<cstdio>
#include<algorithm>
#include<cmath>
#include<bitset>
#pragma warning(disable:4996)
using namespace std;
unsigned HPrime[20000] = { 5,9 }, _PTy, MaxHPrime, * HPrime_end = HPrime + sizeof(HPrime) / sizeof(HPrime[0]);
unsigned x, num[1000002], n = 0; bitset<1000002> isHSemiPrime; bitset<250000> isHPrime;
inline void genprime() {
	unsigned a = 13, t; bool flag = true;
	for (unsigned* i = HPrime + 2; i != HPrime_end;) {
		t = sqrt(a); flag = true;
		for (unsigned* j = HPrime; *j <= t; ++j)
			if (a % *j == 0) { flag = false; break; }
		if (flag) { *i = a, ++i, isHPrime[a] = 1; }
		a += 4;
	}
	MaxHPrime = *(HPrime_end - 1);
}
int main() {
	genprime(); isHPrime[5] = isHPrime[9] = 1;
	for (unsigned i = 5; i <= MaxHPrime; i += 4) {
		if (isHPrime[i] == 0)continue;
		for (unsigned j = 5; j <= MaxHPrime; j += 4) {
			if (isHPrime[j] == 0)continue;
			static unsigned k; k = i * j; if (k > 1000001)break;
			isHSemiPrime[k] = 1;
		}
	}
	for (unsigned i = 5; i <= 1000001; i += 4) {
		if (isHSemiPrime[i] == 1)++n;
		num[i] = n;
	}
	for (;;) {
		scanf("%u", &x); if (x == 0)return 0;
		printf("%u %u\n", x, num[x]);
	}
}