prime

使用加法,徹底分解數字

例如5 = 1 + 1 + 1 + 1 + 1。

不可分解的單元是1。不稀奇。

使用乘法,徹底分解數字

加法換成乘法之後,事情變化多端,例如5 = 5,例如6 = 2 × 3。可以發現,2、3、5、7、11等,是不可分解的單元。至於4、6、8、9、10等,可以再分解。

不可分解的單元叫做「質數」!非常稀奇!

接下來要介紹的演算法有:從小到大列出質數(建立質數表)、判斷一個數是不是質數(質數測試)、使用乘法湊得給定的數(質因數分解)。

prime generation: sieve of Eratosthenes

prime generation

製作質數表。

  2   3   5   7  11  13  17  19  23  29
 31  37  41  43  47  53  59  61  67  71
 73  79  83  89  97 101 103 107 109 113
127 131 137 139 149 151 157 163 167 173
179 181 191 193 197 199 211 223 227 229
233 239 241 251 257 263 269 271 277 281
......

sieve of Eratosthenes

這是一個製作質數表的演算法。簡稱「篩法」。

列出所有正整數。從2開始,刪掉2的倍數。找下一個未被刪掉的數字,找到3,刪掉3的倍數。找下一個未被刪掉的數字,找到5,刪掉5的倍數。如此不斷下去,就能刪掉所有合數,找到所有質數。

質數有無限多個,證明省略。我們無法找到所有質數,通常是預先訂立一個範圍,只找到範圍內的所有質數。

欲刪掉質數i的倍數之時,早已刪掉1倍到i-1倍了,直接從i倍開始。

欲刪掉質數i的倍數之時,早已刪掉「小於i的質數、其倍數」倍了,直接刪掉「大於等於i的質數、其倍數」倍。

乍看下程式碼增多而變慢,實際上cache miss減少而變快。

一個合數x,必定有一個小於等於sqrt(x)的質因數。所有小於等於sqrt(x)的質數,刪掉這些質數的倍數,就能刪掉所有合數了。

顛倒true和false,節省初始化時間。

製作質數表。篩法結束之後,掃描一次陣列即可。

UVa 406 516 524 543 10140 10311 10394 11408

使用bitset來取代bool陣列

一個int有32個位元,可以當作32個欄位來使用,節省記憶體空間,減少cache miss。

不處理2的倍數

不處理2的倍數,節省一半記憶體,增進一點速度。

令陣列第0格代表數字1、陣列第1格代表數字3、陣列第2格代表數字5、……,以此類推。

時間複雜度

考慮內層迴圈索引值j一共有多少種:(N/2 - 2) + (N/3 - 3) + (N/5 - 5) + ... + (N/sqrtN - sqrtN) = O(NloglogN)。

1/1 + 1/2 + 1/3 + ... + 1/N     = O(logN)
1/2 + 1/3 + 1/5 + ... + 1/N     = O(loglogN)
1/2 + 1/3 + 1/5 + ... + 1/sqrtN = O(loglogsqrtN) = O(loglogN)
https://en.wikipedia.org/wiki/Divergence_of_the_sum_of_the_reciprocals_of_the_primes
1/1 + 1/2 + 1/3 + ... + 1/n - ln(n)     趨近於Euler–Mascheroni constant
1/2 + 1/3 + 1/5 + ... + 1/n - ln(ln(n)) 趨近於Meissel–Mertens constant

prime generation: segmented sieve of Eratosthenes

segmented sieve of Eratosthenes

篩法當中,以不大於sqrtN的質數,篩選大於sqrtN的數。

篩法需要大量記憶體空間。為了節省記憶體空間,於是有了分段處理的手法:一、首先求得不大於sqrtN的質數。二、將記憶體切成小段,每段長度sqrtN。分別篩選每一段。

空間複雜度從O(N)變成O(sqrtN),大幅減少cache miss。

也許可以作為平行演算法的經典範例。

prime generation: linear sieve algorithm

一次方時間的篩法

一邊製作質數表,一邊刪掉每個數的質數倍。每個合數只遇到最小質因數(次小因數),每個合數只讀取一次。時間複雜度O(N)。

primality test

primality test

「質數測試」。測試一個數字是否為質數。

質數測試屬於P問題。然而以下介紹的都是指數時間演算法,甚至是不保證結果正確的演算法。

多項式時間、保證結果正確的演算法,請逕行上網搜尋「AKS primality test」;但由於方法繁雜,實務上沒人使用。

質數測試,可以直接以篩法建立質數表,再來判斷質數。然而建立質數表需要大量記憶體,因此又發明了其他方法。

primality test: divisibility primality test

divisibility primality test

整除性測試法即是試除法。

質數無法以乘法分解,質數沒有任何因數(除了1和自身)。窮舉所有可能的因數,一一試除。時間複雜度O(n)。

因數成雙成對(平方根跟自己一對)。檢查了小於等於平方根的因數,就不必再檢查大於平方根的因數。時間複雜度O(sqrt(n))。

2n+1 method

因數可以分解成質數。取質數試除,效果等同於取因數試除。反過來說,大可不必取合數試除。

2的更大倍數是合數,不必試除,節省一半時間。

6n±1 method

2和3的最小公倍數是6。所有數字分為6n、6n+1、6n+2、6n+2、6n+3、6n+4、6n+5六種(n是倍率)。除以六餘零的數字為6n,除以六餘一的數字為6n+1,以此類推。

6n、6n+2、6n+3、6n+4是2或3的倍數,不是質數。只需要拿6n+1和6n+5試除。

從6n+1,跳到6n+5,跳到下一個6n+1,再跳到6n+5,差值是4 2 4 2 4 2...,不斷跳四跳二。

6n+5改寫成6n-1,意義不變。6n+1與6n-1簡寫成6n±1。

2和3預先處理。從5開始試除,下一個是5+2,再下一個是5+2+4,再下一個是5+2+4+2,以此類推。

實作程式碼,用加法加二加四,不用乘法計算6n-1、6n+1,效率較佳。

wheel factorization

6n±1法進行推廣,使用頭幾個質數。

取2 3 5的最小公倍數30,建立30n+k法。取2 3 5 7的最小公倍數210,建立210n+k法。取所有質數,等同篩法。

取越多質數,即可剔除越多合數,事前節省時間。

然而,取越多質數,也會建立越多規則,事中浪費時間。

這是取捨兩難的情況,兩者之間必須取得平衡。一般來說,大家只取2 3或者2 3 5。

更有甚者,當你還在煩惱要取多少質數,普通的試除法早就已經執行完畢了。更有甚者,當你還在煩惱要如何實作試除法,網路上已經有現成的質數表了。所謂的取得平衡,端看你的視野高度。

primality test: Miller–Rabin primality test

square root primality test

餘數系統的平方根有個特性:

以質數n為模數,0到n-1之間,只有1與n-1,平方之後等於1。
以質數n為模數,1的平方根一定是1與n-1(即±1)。

平方根測試法,運用了此特性:

以n為模數,當2到n-2的每一個數字,平方之後皆不等於1,就推定n是質數。

此演算法的結果不一定正確。通過測試的數字,可能是合數或質數;無法通過測試的數字,一定是合數。

以合數n為模數,0到n-1之間,除了1與n-1,還可能有其它數字,平方之後等於1。
以合數n為模數,1的平方根可能不只是1與n-1(不只±1)。

有些合數會被誤判成質數,例如22會被誤判成質數,2²、…、20²模22之後皆不等於1。

Fermat primality test

費瑪小定理:

若n是質數,a小於n,則aⁿ ≡ a (mod n)。
若n是質數,a小於n,a不是零,則aⁿ⁻¹ ≡ 1 (mod n)。

費瑪質數測試法,運用了費瑪小定理:

n是質數,費瑪小定理一定成立:aⁿ⁻¹ % n = 1一定成立。
n是合數,費瑪小定理可能成立:aⁿ⁻¹ % n = 1可能成立。

當aⁿ⁻¹ % n = 1成立,就推定n是質數。

此演算法的結果不一定正確。通過測試的數字,可能是合數或質數;無法通過測試的數字,一定是合數。

使用各式各樣的a來實施測試,那麼判定結果就更準確。

Miller–Rabin primality test

結合square root primality test與Fermat primality test。

此演算法的結果不一定正確。通過測試的數字,可能是合數或質數;無法通過測試的數字,一定是合數。

一、選定一個底數a,大於1、小於n,用來進行費瑪測試。
  待測數字n,會是費瑪測試的模數。
二、令 n-1 = u × 2ᵗ,其中t盡量大(u為奇數)。
三、觀察a^u這個數字:
  若等於±1,
  則表示步驟四,所有數字都是1,永不出現平方根測試的反例。
  也表示步驟五,最終將通過費瑪測試。
  推定n是質數。
四、依序觀察a^(u × 2¹)、a^(u × 2²)、…、a^(u × 2ᵗ⁻¹)這些數字,
  每個數字正好是前一個數字的平方:
 甲、一旦發現有個數字的平方等於+1,
   則表示無法通過平方根測試。
   (但是步驟五,最終將通過費瑪測試。)
   確定n是合數。
 乙、一旦發現有個數字的平方等於-1,
   則表示接下來的數字都是1,永不出現平方根測試的反例。
   也表示步驟五,最終將通過費瑪測試。
   推定n是質數。
五、觀察a^(u × 2ᵗ)這個數字:
 甲、若等於+1,表示通過費瑪測試,推定n是質數。
 乙、若不等於+1,表示無法通過費瑪測試,確定n是合數。
 回、也就是說,a^(u × 2ᵗ⁻¹)必須等於±1,平方之後才會等於+1。
   步驟五才能通過費瑪測試。
 回、由於步驟四已經檢查過a^(u × 2ᵗ⁻¹)是否為±1,
   所以步驟五可以省略,直接確定n是合數。

步驟二:不斷分解因數2,最多log(n)步,需時O(log(n))。

步驟三:餘數次方,divide-and-conquer method,需時O(log(n))。

步驟四:根據步驟二,最多log(n)-1個數字,需時O(log(n))。

總時間複雜度O(log(n))。

strong probable-prime base(SPRP)

以過濾合數的觀點來看,多取幾個相異的底數a實施測試,判定結果就更準確。

事實上已經有熱心人士,找出特定的底數組合,仔細檢查了各種數字的判定結果,保證判定結果正確。例如底數組合{2, 7, 61}就能正確判斷2³²以下的數字是不是質數。

UVa 10956

primality test: Baillie–PSW primality test

Lucas primality test

Baillie–PSW primality test

目前仍然可以正確判斷質數,尚未找到反例。