第一章 素?cái)?shù)判定法現(xiàn)狀
現(xiàn)在，確定性素?cái)?shù)判定法已經(jīng)有很多種，常用的有試除法、威廉斯方法、艾德利曼和魯梅利法。它們的適用范圍各不相同，威廉斯方法比較適合10^20到10^50之間的數(shù)，艾德利曼和魯梅利法適合大于10^50的數(shù)，對于32位機(jī)器數(shù)，由于都小于10^10，所以一般都用試除法來判定。
也許有人會(huì)問：“你為什么沒有提馬寧德拉.阿格拉瓦法呢？不是有人說它是目前最快的素?cái)?shù)判定法嗎？” 其實(shí)這是一個(gè)很大的誤解，阿格拉瓦法雖然是log(n)的多項(xiàng)式級算法，但目前只有理論上的意義，根本無法實(shí)用，因?yàn)樗臅r(shí)間復(fù)雜度是O(log(n)^12)，這個(gè)多項(xiàng)式的次數(shù)太高了。就拿最慢的試除法跟它來比吧，試除法的時(shí)間復(fù)雜度為O(n^(1/2)*log(n)^2)，當(dāng)n = 16時(shí)，log(n)^12 = 16777216，而n^(1/2)*log(n)^2 = 64，你看相差有多么大！如果要讓兩者速度相當(dāng)，即log(n)^12 = n^(1/2)*log(n)^2，得出n = 10^43.1214，此時(shí)需要進(jìn)行的運(yùn)算次數(shù)為log(n)^12 = 10^25.873（注意：本文中log（）函數(shù)缺省以2為底），這樣的運(yùn)算次數(shù)在一臺(tái)主頻3GHz的計(jì)算機(jī)上運(yùn)行也要10^8.89707年才能運(yùn)行完，看來我們這輩子是別指望看到阿格拉瓦法比試除法快的這一天啦！
除了這些確定性素?cái)?shù)判定法外，還有基于概率的非確定性素?cái)?shù)判定法，最常用的就是米勒-拉賓法。
對于32位機(jī)器數(shù)（四則運(yùn)算均為常數(shù)時(shí)間完成），試除法的時(shí)間復(fù)雜度是O(n^(1/2))，而米勒-拉賓法的時(shí)間復(fù)雜度只有O(log(n))。所以后者要比前者快得多，但是由于米勒-拉賓法的非確定性，往往我們在需要確定解時(shí)仍然要依靠速度較慢的試除法。那是否可以通過擴(kuò)展米勒-拉賓法，來找到一種更快的確定性素?cái)?shù)判定法呢？結(jié)論是肯定的，本文就帶你一起尋找這樣一種方法。
第二章 2-偽素?cái)?shù)表的生成
既然要擴(kuò)展米勒-拉賓法，那首先我們應(yīng)該知道為什么米勒-拉賓法是個(gè)非確定性素?cái)?shù)判定法？答案很簡單，由于偽素?cái)?shù)的存在。由于米勒-拉賓法使用費(fèi)爾馬小定理的逆命題進(jìn)行判斷，而該逆命題對極少數(shù)合數(shù)并不成立，從而產(chǎn)生誤判，這些使費(fèi)爾馬小定理的逆命題不成立的合數(shù)就是偽素?cái)?shù)。為了研究偽素?cái)?shù)，我們首先需要生成偽素?cái)?shù)表，原理很簡單，就是先用篩法得出一定范圍內(nèi)的所有素?cái)?shù)，然后逐一判定該范圍內(nèi)所有合數(shù)是否使以2為底數(shù)的費(fèi)爾馬小定理的逆命題不成立，從而得出該范圍內(nèi)的2-偽素?cái)?shù)表。我的程序運(yùn)行了100分鐘，得出了32位機(jī)器數(shù)范圍內(nèi)的2-偽素?cái)?shù)表，如下：
341
561
645
1105
1387
1729
1905
2047
2465
2701
...
...
...
4286813749
4288664869
4289470021
4289641621
4289884201
4289906089
4293088801
4293329041
4294868509
4294901761
（共10403個(gè)，由于篇幅所限，中間部分省略。）
第三章 尋找2-偽素?cái)?shù)的最小可判定底數(shù)
對于2-偽素?cái)?shù)表的每一個(gè)偽素?cái)?shù)，尋找最小的可以判定它們是合數(shù)的底數(shù)，我把這個(gè)底數(shù)稱之為最小可判定底數(shù)。特別地，對于絕對偽素?cái)?shù)，它的最小質(zhì)因子即是它的最小可判定底數(shù)。由于已經(jīng)證明了絕對偽素?cái)?shù)至少有三個(gè)質(zhì)因子，所以這個(gè)最小質(zhì)因子一定不大于n^(1/3)。下面就是我找到的最小可判定底數(shù)列表：
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
通過統(tǒng)計(jì)這個(gè)列表，我發(fā)現(xiàn)了一個(gè)規(guī)律，那就是所有的最小可判定底數(shù)都不大于n^(1/3)，由前述可知，對于絕對偽素?cái)?shù)，這個(gè)結(jié)論顯然成立。而對于非絕對偽素?cái)?shù)，雖然直觀上覺得它應(yīng)該比絕對偽素?cái)?shù)好判定出來，但是我無法證明出它的最小可判定底數(shù)都不大于n^(1/3)。不過沒關(guān)系，這個(gè)問題就作為一個(gè)猜想留給數(shù)學(xué)家來解決吧，更重要的是我已經(jīng)通過實(shí)驗(yàn)證明了在32位機(jī)器數(shù)范圍內(nèi)這個(gè)結(jié)論成立。
我們還有沒有更好的方法來進(jìn)一步減小最小可判定底數(shù)的范圍呢？有的！我們可以在計(jì)算平方數(shù)時(shí)進(jìn)行二次檢測，下面是進(jìn)行了二次檢測后重新計(jì)算的最小可判定底數(shù)列表：
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
很顯然，二次檢測是有效果的，經(jīng)過統(tǒng)計(jì)，我發(fā)現(xiàn)了新的規(guī)律，那就是經(jīng)過二次檢測后所有的最小可判定底數(shù)都不大于n^(1/6)，真的是開了一個(gè)平方呀，哈哈！這個(gè)結(jié)論的數(shù)學(xué)證明仍然作為一個(gè)猜想留給數(shù)學(xué)家們吧。我把這兩個(gè)猜想叫做費(fèi)爾馬小定理可判定上界猜想。而我已經(jīng)完成了對32位機(jī)器數(shù)范圍內(nèi)的證明。
通過上面總結(jié)的規(guī)律，我們已經(jīng)可以設(shè)計(jì)出一個(gè)對32位機(jī)器數(shù)進(jìn)行素?cái)?shù)判定的 O(n^(1/6)*log(n)) 的確定性方法。但是這還不夠，我們還可以優(yōu)化，因?yàn)榇藭r(shí)的最小可判定底數(shù)列表去重后只剩下了5個(gè)數(shù)（都是素?cái)?shù)）：{2,3,5,7,11}。天哪，就是前5個(gè)素?cái)?shù)，這也太容易記憶了吧。
不過在實(shí)現(xiàn)算法時(shí)，需要注意這些結(jié)論都是在2-偽素?cái)?shù)表基礎(chǔ)上得來的，也就是說不管如何對2的判定步驟必不可少，即使當(dāng)2>n^(1/6)時(shí)。
還有一些優(yōu)化可以使用，經(jīng)過實(shí)驗(yàn)，當(dāng)n>=7^6時(shí)，可以不進(jìn)行n^(1/6)上界限制，而固定地用{2,5,7,11}去判定，也是100%正確的。這樣就可以把判定次數(shù)降為4次以下，而每次判定只需要進(jìn)行4log(n)次乘除法（把取余運(yùn)算也看作除法），所以總的計(jì)算次數(shù)不會(huì)超過16log(n)。經(jīng)過實(shí)驗(yàn)，最大的計(jì)算次數(shù)在n=4294967291時(shí)出現(xiàn)，為496次。
自己寫了個(gè)隨機(jī)化素?cái)?shù)判定算法：
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 //統(tǒng)計(jì)二進(jìn)制位數(shù)
16 while (power1 > 0)
17 {
18 power1 >>= 1;
19 bits++;
20 }
21 //從高位到低位依次處理power的二進(jìn)制位
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 //費(fèi)爾馬小定理逆命題判定
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越大，穩(wěn)定性越好，但效率會(huì)低點(diǎn)。原作者有更好的判定算法，對作者表示感謝。。
posted on 2008-01-26 09:37
小果子 閱讀(1202)
評論(0) 編輯 收藏 引用 所屬分類:
學(xué)習(xí)筆記