查看文章 |
现在，确定性素数判定法已经有很多种，常用的有试除法、威廉斯方法、艾德利曼和鲁梅利法。它们的适用范围各不相同，威廉斯方法比较适合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次。 自己写了个随机化素数判定算法： 1 #include <iostream>
s越大，稳定性越好，但效率会低点。原作者有更好的判定算法，对作者表示感谢。。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 } |
