跳轉到

2.1 質數測試與篩法#

概念說明#

質數是大於 1 的整數中，只能被 1 與自身整除者。質數是整數的「原子」：根據算術基本定理，每個大於 1 的整數都能唯一分解成質數的乘積。因此「判斷一個數是不是質數」與「找出某範圍內的所有質數」是數論演算法的兩大基本功。

最直覺的質數測試是試除法：檢查 n 是否有 2 到 √n 之間的因數。若 n 很大（例如幾百位數的密碼學質數），試除法完全不可行，此時需要米勒–拉賓（Miller–Rabin）概率測試：它利用費馬小定理與平方根技巧，每次隨機測試出錯機率至多 1/4，重複數十次後出錯機率小於 2⁻⁶⁰，實務上被當作確定性演算法使用。而要「列出」n 以內所有質數，最快的方法是埃拉托斯特尼篩法（Sieve of Eratosthenes）。

演算法步驟#

埃拉托斯特尼篩法：

  1. 建立布林陣列 is_prime[0..n]，初始全為 True（0、1 除外）。
  2. 從 p = 2 開始，若 is_prime[p] 為 True，則把 p², p²+p, p²+2p, … 全標記為 False（它們都有因數 p）。
  3. p 遞增到 √n 為止；最後仍為 True 的就是所有質數。

米勒–拉賓測試（對奇數 n）：

  1. 寫 n − 1 = 2^s · d，其中 d 為奇數。
  2. 隨機取 a ∈ [2, n−2]，計算 x = a^d mod n；若 x 為 1 或 n−1 則此輪通過。
  3. 否則連續平方至多 s−1 次，若某次得到 n−1 則通過；都沒有則 n 必為合數。
  4. 重複多輪不同的 a，全部通過則「很可能」是質數。

Python 程式碼範例#

import random

def sieve(n):
    """回傳 0..n 中所有質數"""
    is_prime = [True] * (n + 1)
    is_prime[0] = is_prime[1] = False
    p = 2
    while p * p <= n:                 # 只需篩到 √n
        if is_prime[p]:
            for m in range(p * p, n + 1, p):
                is_prime[m] = False   # 從 p² 開始標記
        p += 1
    return [i for i, flag in enumerate(is_prime) if flag]

def is_prime_miller_rabin(n, rounds=20):
    """米勒-拉賓概率質數測試"""
    if n < 2:
        return False
    for p in (2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37):
        if n % p == 0:
            return n == p
    d, s = n - 1, 0
    while d % 2 == 0:
        d //= 2
        s += 1
    for _ in range(rounds):
        a = random.randrange(2, n - 1)
        x = pow(a, d, n)              # 模冪（見 2.3 節）
        if x in (1, n - 1):
            continue
        for _ in range(s - 1):
            x = x * x % n
            if x == n - 1:
                break
        else:
            return False              # 一定是合數
    return True

print(sieve(50))
# [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
print(is_prime_miller_rabin(2**61 - 1))       # True（梅森質數）
print(is_prime_miller_rabin(2**61 - 1 + 2))   # False
print(is_prime_miller_rabin(10**18 + 9))      # 大數也能瞬間判斷

複雜度分析#

試除法是 O(√n) 次除法，篩法是 O(n log log n) 時間、O(n) 空間——log log n 來自「對每個質數 p，標記次數為 n/p，對所有 p 求和得 n·Σ(1/p) ≈ n·ln ln n」。從 p² 開始標記與只篩到 √n 是兩個重要常數最佳化。米勒–拉賓每輪需要 O(log n) 次模乘法，20 輪總計約 20·log n 次，即使 n 是 200 位數也只需毫秒級（Python 內建 pow 用 C 實作大數模乘）。Python 3.8+ 其實也提供 sympy.isprime，內部就是確定性的米勒–拉賓變體，可用來驗證自己的實作。

練習#

  1. 修改 sieve，改用「奇數專用陣列」把記憶體砍半，並比較 n = 10⁷ 時兩者的執行時間。
  2. 實作「線性篩」（歐拉篩）：每個合數只被其最小質因數標記一次，時間嚴格 O(n)。提示：維護已找到的質數列表。
  3. 用 is_prime_miller_rabin 找出 10¹⁸ 之後的第一個質數，並與 sympy.nextprime(10**18) 比對。