一种适合 32 位机器数的确定性素数判定法
作者：
王浩（great_wh@163.com ）
此论文是为发布个人研究成果与大家分享交流而提交，本人教育背景如下：
计算机应用专业 / 管理工程专业双学士
天津大学
于 1995 年到 1999 年在校就读
日期________________2006-11-28___________________
摘要
一种适合 32 位机器数的确定性素数判定法
作者： 王浩
本文通过对 米勒- 拉宾非确定性素数判定法如何转化为确定性素数判定法的研究，发现了与之相关的伪素数的一些性质，引入了伪素数的最小可判定底数的概念，并总结出了一些规律。通过这些规律找出了 一种特别适合 32 位机器数的确定性素数判定法，该方法对于 32 位机器数进行素数判定最多只需要进行 16log(n) 次乘 / 除法。该方法具有实现简单、速度快的优点，非常具有推广价值。
本文中总结出的一些规律如果能够得到证明和推广，则有可能彻底解决把米勒- 拉宾非确定性素数判定法转化为确定性素数判定法的问题，从而对素数判定理论和实践产生一定的促进作用。
本文共有五章。分述如下：
第一章：讲述素数判定法的现状，列举了目前常用的一些素数判定法及其适用范围。
第二章：讲解伪素数表生成过程。
第三章：分析伪素数表，引入了伪素数的最小可判定底数的概念，并且总结出了一些规律。根据这些规律，找出了 一种特别适合 32 位机器数的确定性素数判定法，并且进行了多种优化，给出了时间复杂度分析。
第四章：算法的C++ 语言实现和解释说明。
第五章：算法的可推广性分析和未来发展展望。
目录
词汇表
素数判定法 ：判定一个自然数是否素数的方法。
确定性素数判定法 ：一个素数判定法判定某个自然数为素数的充要条件是该自然数确实是素数，该判定法就是确定性素数判定法。即该判定法不存在误判的可能性。
32 位机器数 ：在计算机上用32 个二进制位表示的无符号整数。
64 位机器数 ：在计算机上用64 个二进制位表示的无符号整数。
现在，确定性素数判定法已经有很多种，常用的有试除法、威廉斯方法、艾德利曼和鲁梅利法。它们的适用范围各不相同，威廉斯方法比较适合10^20 到10^50 之间的数，艾德利曼和鲁梅利法适合大于10^50 的数，对于32 位机器数，由于都小于10^10 ，所以一般都用试除法来判定。
也许有人会问：“你为什么没有提马宁德拉. 阿格拉瓦法呢？不是有人说它是目前最快的素数判定法吗？” 其实这是一个很大的误解，阿格拉瓦法虽然是log(n) 的多项式级算法，但目前只有理论上的意义，根本无法实用，因为它的时间复杂度是O(log(n)^12) ，这个多项式的次数太高了。就拿最慢的试除法跟它来比吧，试除法的时间复杂度为O(n^(1/2)*log(n)^2) ，当n = 16 时，log(n)^12 = 16777216 ，而n^(1/2)*log(n)^2 = 64 ，你看相差有多么大！如果要让两者速度相当，即log(n)^12 = n^(1/2)*log(n)^2 ，得出n = 10^43.1214 ，此时需要进行的运算次数为log(n)^12 = 10^25.873 （注意：本文中log （）函数缺省以2 为底），这样的运算次数在一台主频3GHz 的计算机上运行也要10^8.89707 年才能运行完，看来我们这辈子是别指望看到阿格拉瓦法比试除法快的这一天啦！
除了这些确定性素数判定法外，还有基于概率的非确定性素数判定法，最常用的就是米勒- 拉宾法。
对于32 位机器数（四则运算均为常数时间完成），试除法的时间复杂度是O(n^(1/2)) ，而米勒- 拉宾法的时间复杂度只有O(log(n)) 。所以后者要比前者快得多，但是由于米勒- 拉宾法的非确定性，往往我们在需要确定解时仍然要依靠速度较慢的试除法。那是否可以通过扩展米勒- 拉宾法，来找到一种更快的确定性素数判定法呢？结论是肯定的，本文就带你一起寻找这样一种方法。
既然要扩展米勒- 拉宾法，那首先我们应该知道为什么米勒- 拉宾法是个非确定性素数判定法？答案很简单，由于伪素数的存在。由于米勒- 拉宾法使用费尔马小定理的逆命题进行判断，而该逆命题对极少数合数并不成立，从而产生误判，这些使费尔马小定理的逆命题不成立的合数就是伪素数。为了研究伪素数，我们首先需要生成伪素数表，原理很简单，就是先用筛法得出一定范围内的所有素数，然后逐一判定该范围内所有合数是否使以2 为底数的费尔马小定理的逆命题不成立，从而得出该范围内的2- 伪素数表。我的程序运行了100 分钟，得出了32 位机器数范围内的2- 伪素数表，如下：
341
561
645
1105
1387
1729
1905
2047
2465
2701
...
...
...
4286813749
4288664869
4289470021
4289641621
4289884201
4289906089
4293088801
4293329041
4294868509
4294901761
（共 10403 个，由于篇幅所限，中间部分省略。）
对于 2- 伪素数表的每一个伪素数，寻找最小的可以判定它们是合数的底数，我把这个底数称之为最小可判定底数。特别地，对于绝对伪素数，它的最小质因子即是它的最小可判定底数。由于已经证明了绝对伪素数至少有三个质因子，所以这个最小质因子一定不大于 n^(1/3) 。下面就是我找到的最小可判定底数列表：
341 3
561 3
645 3
1105 5
1387 3
1729 7
1905 3
2047 3
2465 5
2701 5
...
...
...
4286813749 3
4288664869 3
4289470021 5
4289641621 3
4289884201 3
4289906089 3
4293088801 3
4293329041 3
4294868509 7
4294901761 3
通过统计这个列表，我发现了一个规律，那就是所有的最小可判定底数都不大于 n^(1/3) ，由前述可知，对于绝对伪素数，这个结论显然成立。而对于非绝对伪素数，虽然直观上觉得它应该比绝对伪素数好判定出来，但是我无法证明出它的最小可判定底数都不大于 n^(1/3) 。不过没关系，这个问题就作为一个猜想留给数学家来解决吧，更重要的是我已经通过实验证明了在 32 位机器数范围内这个结论成立。
我们还有没有更好的方法来进一步减小最小可判定底数的范围呢？有的！我们可以在计算平方数时进行二次检测，下面是进行了二次检测后重新计算的最小可判定底数列表：
341 2
561 2
645 2
1105 2
1387 2
1729 2
1905 2
2047 3
2465 2
2701 2
...
...
...
4286813749 2
4288664869 2
4289470021 2
4289641621 2
4289884201 2
4289906089 2
4293088801 2
4293329041 2
4294868509 2
4294901761 3
很显然，二次检测是有效果的，经过统计，我发现了新的规律，那就是经过二次检测后所有的最小可判定底数都不大于 n^(1/6) ，真的是开了一个平方呀，哈哈！这个结论的数学证明仍然作为一个猜想留给数学家们吧。我把这两个猜想叫做费尔马小定理可判定上界猜想。而我已经完成了对 32 位机器数范围内的证明。
通过上面总结的规律，我们已经可以设计出一个对 32 位机器数进行素数判定的 O(n^(1/6)*log(n)) 的确定性方法。但是这还不够，我们还可以优化，因为此时的最小可判定底数列表去重后只剩下了 5 个数（都是素数）： {2,3,5,7,11} 。天哪，就是前 5 个素数，这也太容易记忆了吧。
不过在实现算法时，需要注意这些结论都是在 2- 伪素数表基础上得来的，也就是说不管如何对 2 的判定步骤必不可少，即使当 2>n^(1/6) 时。
还有一些优化可以使用，经过实验，当 n>=7^6 时，可以不进行 n^(1/6) 上界限制，而固定地用 {2,5,7,11} 去判定，也是 100% 正确的。这样就可以把判定次数降为 4 次以下，而每次判定只需要进行 4log(n) 次乘除法（把取余运算也看作除法），所以总的计算次数不会超过 16log(n) 。经过实验，最大的计算次数在 n=4294967291 时出现，为 496 次。
算法实现如下：（使用 C++ 语言）
#include <iostream>
#include <math.h>
using namespace std;
// 定义跨平台的 64 位机器数类型
#ifndef _WIN32
typedef unsigned long long longlong_t;
#else
typedef unsigned __int64 longlong_t;
#endif
// 使用费尔马小定理和二次检测针对一个底数进行判定
bool IsLikePrime(longlong_t n, longlong_t base)
{
longlong_t power = n-1;
longlong_t result = 1;
longlong_t x = result;
longlong_t bits = 0;
longlong_t power1 = power;
// 统计二进制位数
while (power1 > 0)
{
power1 >>= 1;
bits++;
}
// 从高位到低位依次处理 power 的二进制位
while(bits > 0)
{
bits--;
result = (x*x)%n;
// 二次检测
if (result == 1 && x != 1 && x != n-1)
{
return false;
}
if ((power&((longlong_t)1<<bits)) != 0)
{
result = (result*base)%n;
}
x = result;
}
// 费尔马小定理逆命题判定
return result == 1;
}
// 前 5 个素数
const int primes[]={2,3,5,7,11};
// 前 5 个素数的 6 次方，由后面的 init 对象初始化
int primes_six[sizeof(primes)/sizeof(primes[0])];
// 静态初始化类
class CInit
{
public:
CInit()
{
int num = sizeof(primes)/sizeof(primes[0]);
for (int i = 0; i < num; i++)
{
primes_six[i] = primes[i]*primes[i]*primes[i];
primes_six[i] *= primes_six[i];
}
}
}init;
// 王浩素数判定函数
bool JudgePrime(longlong_t n)
{
if (n < 2)
return false;
if (n == 2)
return true;
int num = sizeof(primes)/sizeof(int);
bool bIsLarge = (n >= primes_six[3]);
for (int i = 0; i < num; i++)
{
if (bIsLarge)
{
// 当 n >= 7^6 时，不进行上界判断，固定地用 {2,5,7,11} 做判定。
if (primes[i] == 3)
continue;
}
else
{
// 当 n < 7^6 时，进行上界判断，但是 2 例外。
if (primes[i] != 2 && n < primes_six[i])
break;
}
// 做一次子判定
if (!IsLikePrime(n, primes[i]))
return false;
}
// 所有子判定通过，则 n 必为素数！
return true;
}
// 主程序
int main()
{
longlong_t n;
// 对标准输入的每一个数进行素数判定
while (cin >> n)
{
if (JudgePrime(n))
{
// 如果是素数，则输出到标准输出。
cout << n << endl;
}
// 如果是合数，不输出。
}
return 0;
}
程序中已经加了足够的注释，应该不难理解。
需要说明的一点是，虽然我在输入时使用了 longlong_t ，那是为了类型一致性，有效的输入范围仍然是 0 ~ 2^32-1 。
如果前述的费尔马小定理可判定上界猜想可以被证明，那么该算法可以被推广到任意位数的 n ，此时的时间复杂度为 O(n^(1/6)*log(n)^3) 。这样我们就可以完成 米勒- 拉宾非确定性素数判定法向确定性素数判定法的转化，这对于数论理论是一个补充，对于实践中使用米勒- 拉宾素数判定法具有指导意义。
本文所做的研究只是向 米勒- 拉宾非确定性素数判定法的确定化方向迈出了一小步，我相信，在不久的将来，米勒- 拉宾非确定性素数判定法的确定化方向会有更大进展，从而对数论理论和实践产生深远影响。
《计算机算法设计与分析（第 2 版）》，王晓东编著，电子工业出版社， 2004 年 7 月。
扫一扫
194
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
