第一章 素数判定法现状
现在，确定性素数判定法已经有很多种，常用的有试除法、威廉斯方法、艾德利曼和鲁梅利法。它们的适用范围各不相同，威廉斯方法比较适合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为底数的费尔马小定理的逆命题不成立，从而得出该范围内的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-伪素数的最小可判定底数
对于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次。
自己写了个随机化素数判定算法：
1 #include <iostream>
2 #include <cmath>
3 #include <algorithm>
4
5 using namespace std;
6 typedef unsigned __int64 longlong_t;
7
8 bool _IsLikePrime(longlong_t n, longlong_t base)
9 {
10 longlong_t power = n-1;
11 longlong_t result = 1;
12 longlong_t x = result;
13 longlong_t bits = 0;
14 longlong_t power1 = power;
15 //统计二进制位数
16 while (power1 > 0)
17 {
18 power1 >>= 1;
19 bits++;
20 }
21 //从高位到低位依次处理power的二进制位
22 while(bits > 0)
23 {
24 bits--;
25 result = (x*x)%n;
26 //二次检测
27 if(result==1&&x!=1&&x!=n-1)
28 return false;
29
30 if ((power&((longlong_t)1<<bits)) != 0)
31 result = (result*base)%n;
32
33 x = result;
34 }
35 //费尔马小定理逆命题判定
36 return result == 1;
37 }
38 inline bool JudgePrime(int n,int s)
39 {
40 srand(20);
41 for(int i=0;i<s;i++){
42 int a=rand()%(n-1)+1;
43 if(!_IsLikePrime(n,a))
44 return false;
45 }
46 return true;
47 }
48 int main()
49 {
50 int n;
51 while(scanf("%d",&n)!=EOF){
52 int num=0;
53 int cnt=0;
54 for(int i=0;i<n;i++){
55 scanf("%d",&num);
56 if(JudgePrime(num,20))
57 cnt++;
58 }
59 printf("%d\n",cnt);
60 }
61 return 0;
62 }
s越大，稳定性越好，但效率会低点。原作者有更好的判定算法，对作者表示感谢。。
posted on 2008-01-26 09:37
小果子 阅读(1174)
评论(0) 编辑 收藏 引用 所属分类:
学习笔记