素數

素數與合數的定義，見數論基礎。

素數計數函數：小於或等於 \(x\) 的素數的個數，用 \(\pi(x)\) 表示。隨着 \(x\) 的增大，有這樣的近似結果：\(\pi(x) \sim \dfrac{x}{\ln(x)}\)。

素數判定

我們自然地會想到，如何用計算機來判斷一個數是不是素數呢？

實現

暴力做法自然可以枚舉從小到大的每個數看是否能整除

C++Python

bool isPrime(a) {
  if (a < 2) return 0;
  for (int i = 2; i < a; ++i)
    if (a % i == 0) return 0;
  return 1;
}

def isPrime(a):
    if a < 2:
        return False
    for i in range(2, a):
        if a % i == 0:
            return False
    return True

這樣做是十分穩妥了，但是真的有必要每個數都去判斷嗎？

很容易發現這樣一個事實：如果 \(x\) 是 \(a\) 的約數，那麼 \(\frac{a}{x}\) 也是 \(a\) 的約數。

這個結論告訴我們，對於每一對 \((x, \frac{a}{x} )\)，只需要檢驗其中的一個就好了。為了方便起見，我們之考察每一對裏面小的那個數。不難發現，所有這些較小數就是 \([1, \sqrt{a}]\) 這個區間裏的數。

由於 \(1\) 肯定是約數，所以不檢驗它。

C++Python

bool isPrime(a) {
  if (a < 2) return 0;
  for (int i = 2; i * i <= a; ++i)
    if (a % i == 0) return 0;
  return 1;
}

def isPrime(a):
    if a < 2:
        return False
    for i in range(2, int(sqrt(a)) + 1):
        if a % i == 0:
            return False
    return True

素性測試

定義

素性測試（Primality test）是一類在 不對給定數字進行素數分解（prime factorization）的情況下，測試其是否為素數的算法。

素性測試有兩種：

確定性測試：絕對確定一個數是否為素數。常見示例包括 Lucas–Lehmer 測試和橢圓曲線素性證明。
概率性測試：通常比確定性測試快很多，但有可能（儘管概率很小）錯誤地將合數識別為質數（儘管反之則不會）。因此，通過概率素性測試的數字被稱為 可能素數，直到它們的素數可以被確定性地證明。而通過測試但實際上是合數的數字則被稱為 偽素數。有許多特定類型的偽素數，最常見的是費馬偽素數，它們是滿足費馬小定理的合數。概率性測試的常見示例包括 Miller–Rabin 測試。

接下來我們將着重介紹幾個概率性素性測試：

Fermat 素性測試

Fermat 素性檢驗 是最簡單的概率性素性檢驗。

我們可以根據費馬小定理得出一種檢驗素數的思路：

基本思想是不斷地選取在 \([2, n-1]\) 中的基 \(a\)，並檢驗是否每次都有 \(a^{n-1} \equiv 1 \pmod n\)。

實現

C++Python

bool millerRabin(int n) {
  if (n < 3) return n == 2;
  // test_time 為測試次數,建議設為不小於 8
  // 的整數以保證正確率,但也不宜過大,否則會影響效率
  for (int i = 1; i <= test_time; ++i) {
    int a = rand() % (n - 2) + 2;
    if (quickPow(a, n - 1, n) != 1) return 0;
  }
  return 1;
}

def millerRabin(n):
    if n < 3:
        return n == 2
    # test_time 為測試次數,建議設為不小於 8
    # 的整數以保證正確率,但也不宜過大,否則會影響效率
    for i in range(1, test_time + 1):
        a = random.randint(0, 32767) % (n - 2) + 2
        if quickPow(a, n - 1, n) != 1:
            return False
    return True

如果 \(a^{n−1} \equiv 1 \pmod n\) 但 \(n\) 不是素數，則 \(n\) 被稱為以 \(a\) 為底的 偽素數。我們在實踐中觀察到，如果 \(a^{n−1} \equiv 1 \pmod n\)，那麼 \(n\) 通常是素數。但這裏也有個反例：如果 \(n = 341\) 且 \(a = 2\)，即使 \(341 = 11 \cdot 31\) 是合數，有 \(2^{340}\equiv 1 {\pmod {341}}\)。事實上，\(341\) 是最小的偽素數基數。

很遺憾，費馬小定理的逆定理並不成立，換言之，滿足了 \(a^{n-1} \equiv 1 \pmod n\)，\(n\) 也不一定是素數。甚至有些合數 \(n\) 滿足對任意滿足 \(n\nmid a\) 的整數 \(a\) 均有 \(a^{n−1} \equiv 1 \pmod n\)，這樣的數稱為 Carmichael 數。

Carmichael 函數

對正整數 \(n\)，定義 Carmichael 函數（卡邁克爾函數）為對任意滿足 \((a,n)=1\) 的整數 \(a\)，使

\[ a^m\equiv 1\pmod n \]

恆成立的最小正整數 \(m\).

即：

\[ \lambda(n)=\max\{\delta_n(a):(a,n)=1\} \]

Carmichael 函數有如下性質：

（Carmichael 定理）對任意素數 \(p\) 和任意正整數 \(r\)，

\[ \lambda\left(p^r\right)=\begin{cases} \frac{1}{2}\varphi\left(p^r\right), & p=2 \land r\geq 3, \\ \varphi\left(p^r\right), & \text{otherwise}. \end{cases} \]

證明

該定理等價於：

若模 \(n=p^r\) 有原根，則 \(\lambda(n)=\varphi(n)\)，否則 \(\lambda(n)=\dfrac{1}{2}\varphi(n)\).

當模 \(p^r\) 有原根時，由原根存在定理可知命題成立。否則 \(p=2\) 且 \(r\geq 3\)，我們有：

\[ \lambda\left(2^r\right)\mid 2^{r-2} \]

又由 \(5^{2^{r-3}}\equiv 1+2^{r-1}\pmod{2^{r-2}}\) 知 \(\lambda\left(2^r\right)>2^{r-3}\)，因此

\[ \lambda\left(p^r\right)=2^{r-2}=\frac{1}{2}\varphi\left(p^r\right) \]

進而有：
1. 對任意正整數 \(n\)，有 \(\lambda(n)\mid \varphi(n)\)
2. 對任意正整數 \(a\)，\(b\)，有 \(a\mid b\implies \lambda(a)\mid \lambda(b)\)
令 \(n\) 的唯一分解式為 \(n=\prod_{i=1}^k p_i^{r_i}\)，則

\[ \lambda(n)=\left[\lambda\left(p_1^{r_1}\right),\lambda\left(p_2^{r_2}\right),\dots,\lambda\left(p_k^{r_k}\right)\right] \]

由中國剩餘定理和 Carmichael 定理易證。

進而有：
1. 對任意正整數 \(a\)，\(b\)，有 \(\lambda([a,b])=[\lambda(a),\lambda(b)]\)

Carmichael 數

對於合數 \(n\)，如果對於所有正整數 \(a\)，\(a\) 和 \(n\) 互素，都有同餘式 \(a^{n-1} \equiv 1 \pmod n\) 成立，則合數 \(n\) 為 Carmichael 數（卡邁克爾數，OEIS:A002997）。

比如 \(561 = 3 \times 11 \times 17\) 就是一個 Carmichael 數，同時也是最小的 Carmichael 數。

我們可以用如下方法判斷合數 \(n\) 是否為 Carmichael 數：

Korselt 判別法¹

合數 \(n\) 是 Carmichael 數當且僅當 \(n\) 無平方因子且對 \(n\) 的任意質因子 \(p\) 均有 \(p-1 \mid n-1\).

上述判別法可簡化為：

Note

合數 \(n\) 是 Carmichael 數當且僅當 \(\lambda(n)\mid n-1\)，其中 \(\lambda(n)\) 為 Carmichael 函數。

Carmichael 數有如下性質：

Carmichael 數無平方因子且至少有 \(3\) 個不同的質因子。
設 \(C(n)\) 為小於 \(n\) 的 Carmichael 數個數，則：
1. （Alford, Granville, Pomerance. 1994²）\(C(n)>n^{2/7}\)
  
  由此可知 Carmichael 數有無限多個。
2. （Erdős. 1956³）\(C(n)<n\exp\left(-c\dfrac{\ln n\ln\ln\ln n}{\ln\ln n}\right)\)，其中 \(c\) 為常數。
  
  由此可知 Carmichael 數的分佈十分稀疏。實際上 \(C(10^9)=646\)，\(C(10^{18})=1~401~644\)⁴.

注意

「若 \(n\) 為 Carmichael 數，則 \(2^n-1\) 也為 Carmichael 數」是錯誤的。

如 \(561=3 \cdot 11 \cdot 17\) 為 Carmichael 數，考慮 \(2^{561}-1\)。

注意到 \(23\cdot 89=2^{11}-1\mid 2^{561}-1\)，由 Korselt 判別法知，若 \(2^{561}-1\) 是 Carmichael 數，則 \(22\) 和 \(88\) 均為 \(2^{561}-2\) 的因子。

而 \(v_2\left(2^{561}-2\right)=1<v_2(88)=3\)，故 \(88\nmid 2^{561}-2\)，因此 \(2^{561}-1\) 不是 Carmichael 數。

Miller–Rabin 素性測試

Miller–Rabin 素性測試（Miller–Rabin primality test）是進階的素數判定方法。它是由 Miller 和 Rabin 二人根據費馬小定理的逆定理（費馬測試）優化得到的。因為和許多類似算法一樣，它是使用偽素數的概率性測試，我們必須使用慢得多的確定性算法來保證素性。然而，實際上沒有已知的數字通過了高級概率性測試（例如 Miller–Rabin）但實際上卻是複合的。因此我們可以放心使用。

在不考慮乘法的複雜度時，對數 \(n\) 進行 \(k\) 輪測試的時間複雜度是 \(O(k \log n)\)。Miller-Rabbin 素性測試常用於對高精度數進行測試，此時時間複雜度是 \(O(k \log^3n)\)，利用 FFT 等技術可以優化到 \(O(k \log^2n \log \log n \log \log \log n)\)。

二次探測定理

如果 \(p\) 是奇素數，則 \(x^2 \equiv 1 \pmod p\) 的解為 \(x \equiv 1 \pmod p\) 或者 \(x \equiv p - 1 \pmod p\)。

要證明該定理，只需將上面的方程移項，再使用平方差公式，得到 \((x+1)(x-1) \equiv 0 \bmod p\)，即可得出上面的結論。

實現

根據卡邁克爾數的性質，可知其一定不是 \(p^e\)。

不妨將費馬小定理和二次探測定理結合起來使用：

將 \(a^{n-1} \equiv 1 \pmod n\) 中的指數 \(n−1\) 分解為 \(n−1=u \times 2^t\)，在每輪測試中對隨機出來的 \(a\) 先求出 \(v = a^{u} \bmod n\)，之後對這個值執行最多 \(t\) 次平方操作，若發現非平凡平方根時即可判斷出其不是素數，否則再使用 Fermat 素性測試判斷。

還有一些實現上的小細節：

對於一輪測試，如果某一時刻 \(a^{u \times 2^s} \equiv n-1 \pmod n\)，則之後的平方操作全都會得到 \(1\)，則可以直接通過本輪測試。
如果找出了一個非平凡平方根 \(a^{u \times 2^s} \not\equiv n-1 \pmod n\)，則之後的平方操作全都會得到 \(1\)。可以選擇直接返回 false，也可以放到 \(t\) 次平方操作後再返回 false。

這樣得到了較正確的 Miller Rabin：（來自 fjzzq2002）

C++Python

bool millerRabin(int n) {
  if (n < 3 || n % 2 == 0) return n == 2;
  int u = n - 1, t = 0;
  while (u % 2 == 0) u /= 2, ++t;
  // test_time 為測試次數，建議設為不小於 8
  // 的整數以保證正確率，但也不宜過大，否則會影響效率
  for (int i = 0; i < test_time; ++i) {
    int a = rand() % (n - 2) + 2, v = quickPow(a, u, n);
    if (v == 1) continue;
    int s;
    for (s = 0; s < t; ++s) {
      if (v == n - 1) break;  // 得到平凡平方根 n-1，通過此輪測試
      v = (long long)v * v % n;
    }
    // 如果找到了非平凡平方根，則會由於無法提前 break; 而運行到 s == t
    // 如果 Fermat 素性測試無法通過，則一直運行到 s == t 前 v 都不會等於 -1
    if (s == t) return 0;
  }
  return 1;
}

def millerRabin(n):
    if n < 3 or n % 2 == 0:
        return n == 2
    u, t = n - 1, 0
    while u % 2 == 0:
        u = u // 2
        t = t + 1
    # test_time 為測試次數,建議設為不小於 8
    # 的整數以保證正確率,但也不宜過大,否則會影響效率
    for i in range(test_time):
        a = random.randint(2, n - 1)
        v = pow(a, u, n)
        if v == 1:
            continue
        s = 0
        while s < t:
            if v == n - 1:
                break
            v = v * v % n
            s = s + 1
        # 如果找到了非平凡平方根，則會由於無法提前 break; 而運行到 s == t
        # 如果 Fermat 素性測試無法通過，則一直運行到 s == t 前 v 都不會等於 -1
        if s == t:
            return False
    return True

另外，假設廣義 Riemann 猜想（generalized Riemann hypothesis, GRH）成立，則對數 \(n\) 最多隻需要測試 \([2, \min\{n-2, \lfloor 2\ln^2 n \rfloor\}]\) 中的全部整數即可確定數 \(n\) 的素性，證明參見注釋 7。

而在 OI 範圍內，通常都是對 \([1, 2^{64})\) 範圍內的數進行素性檢驗。對於 \([1, 2^{32})\) 範圍內的數，選取 \(\{2, 7, 61\}\) 三個數作為基底進行 Miller–Rabin 素性檢驗就可以確定素性；對於 \([1, 2^{64})\) 範圍內的數，選取 \(\{2, 325, 9375, 28178, 450775, 9780504, 1795265022\}\) 七個數作為基底進行 Miller–Rabin 素性檢驗就可以確定素性。參見注釋 8。

也可以選取 \(\{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37\}\)（即前 \(12\) 個素數）檢驗 \([1, 2^{64})\) 範圍內的素數。

注意如果要使用上面的數列中的數 \(a\) 作為基底判斷 \(n\) 的素性：

所有的數都要取一遍，不能只選小於 \(n\) 的；
把 \(a\) 換成 \(a \bmod n\)；
如果 \(a \equiv 0 \pmod n\)，則直接通過該輪測試。

反素數

引入

顧名思義，素數就是因子只有兩個的數，那麼反素數，就是因子最多的數（並且因子個數相同的時候值最小），所以反素數是相對於一個集合來説的。

一種符合直覺的反素數定義是：在一個正整數集合中，因子最多並且值最小的數，就是反素數。

定義

如果某個正整數 \(n\) 滿足如下條件，則稱為是 反素數：任何小於 \(n\) 的正數的約數個數都小於 \(n\) 的約數個數。

注意

注意區分 emirp，它表示的是逐位反轉後是不同素數的素數（如 149 和 941 均為 emirp，101 不是 emirp）。

過程

那麼，如何來求解反素數呢？

首先，既然要求因子數，首先要做的就是素因子分解。把 \(n\) 分解成 \(n=p_{1}^{k_{1}}p_{2}^{k_{2}} \cdots p_{n}^{k_{n}}\) 的形式，其中 \(p\) 是素數，\(k\) 為他的指數。這樣的話總因子個數就是 \((k_1+1) \times (k_2+1) \times (k_3+1) \cdots \times (k_n+1)\)。

但是顯然質因子分解的複雜度是很高的，並且前一個數的結果不能被後面利用。所以要換個方法。

我們來觀察一下反素數的特點。

反素數肯定是從 \(2\) 開始的連續素數的冪次形式的乘積。
數值小的素數的冪次大於等於數值大的素數，即 \(n=p_{1}^{k_{1}}p_{2}^{k_{2}} \cdots p_{n}^{k_{n}}\) 中，有 \(k_1 \geq k_2 \geq k_3 \geq \cdots \geq k_n\)

解釋：

如果不是從 \(2\) 開始的連續素數，那麼如果冪次不變，把素數變成數值更小的素數，那麼此時因子個數不變，但是 \(n\) 的數值變小了。交換到從 \(2\) 開始的連續素數的時候 \(n\) 值最小。
如果數值小的素數的冪次小於數值大的素數的冪，那麼如果把這兩個素數交換位置（冪次不變），那麼所得的 \(n\) 因子數量不變，但是 \(n\) 的值變小。

另外還有兩個問題，

對於給定的 \(n\)，要枚舉到哪一個素數呢？

最極端的情況大不了就是 \(n=p_{1}p_{2} \cdots p_{n}\)，所以只要連續素數連乘到剛好小於等於 \(n\) 就可以的呢。再大了，連全都一次冪，都用不了，當然就是用不到的啦！
我們要枚舉到多少次冪呢？

我們考慮一個極端情況，當我們最小的素數的某個冪次已經比所給的 \(n\)（的最大值）大的話，那麼展開成其他的形式，最大冪次一定小於這個冪次。unsigned long long 的最大值是 \(2^{64} - 1\)，所以可以枚舉到 \(2^{64} - 1\)。

細節有了，那麼我們具體如何具體實現呢？

我們可以把當前走到每一個素數前面的時候列舉成一棵樹的根節點，然後一層層的去找。找到什麼時候停止呢？

當前走到的數字已經大於我們想要的數字了
當前枚舉的因子已經用不到了（和 \(1\) 重複了嘻嘻嘻）
當前因子大於我們想要的因子了
當前因子正好是我們想要的因子（此時判斷是否需要更新最小 \(ans\)）

然後 dfs 裏面不斷一層一層枚舉次數繼續往下迭代可以。

常見題型

求因子數一定的最小數

例題 Codeforces 27E. A number with a given number of divisors

求具有給定除數的最小自然數。請確保答案不超過 \(10^{18}\)。

解題思路

對於這種題，我們只要以因子數為 dfs 的返回條件基準，不斷更新找到的最小值就可以了

參考代碼

#include <stdio.h>
unsigned long long p[16] = {
    2,  3,  5,  7,  11, 13, 17, 19,
    23, 29, 31, 37, 41, 43, 47, 53};  // 根据数据范围可以确定使用的素数最大为53

unsigned long long ans;
unsigned long long n;

// depth: 当前在枚举第几个素数
// temp: 当前因子数量为 num的时候的数值
// num: 当前因子数
// up：上一个素数的幂，这次应该小于等于这个幂次嘛
void dfs(unsigned long long depth, unsigned long long temp,
         unsigned long long num, unsigned long long up) {
  if (num > n || depth >= 16) return;  // 边界条件
  if (num == n && ans > temp) {        // 取最小的ans
    ans = temp;
    return;
  }
  for (int i = 1; i <= up; i++) {
    if (temp * p[depth] > ans)
      break;  // 剪枝：如果加一个这个乘数的结果比ans要大，则必不是最佳方案
    dfs(depth + 1, temp = temp * p[depth], num * (i + 1),
        i);  // 取一个该乘数，进行对下一个乘数的搜索
  }
}

int main() {
  scanf("%llu", &n);
  ans = ~(unsigned long long)0;
  dfs(0, 1, 1, 64);
  printf("%llu\n", ans);
  return 0;
}

求 n 以內因子數最多的數

例題 ZOJ - More Divisors

大家都知道我們使用十進制記數法，即記數的基數是 \(10\)。歷史學家説這是因為人有十個手指，也許他們是對的。然而，這通常不是很方便，十隻有四個除數——\(1\)、\(2\)、\(5\) 和 \(10\)。因此，像 \(\frac{1}{3}\)、\(\frac{1}{4}\) 或 \(\frac{1}{6}\) 這樣的分數不便於用十進制表示。從這個意義上説，以 \(12\)、\(24\) 甚至 \(60\) 為底會方便得多。主要原因是這些數字的除數要大得多——分別是 \(6\)、\(8\) 和 \(12\)。請回答：除數最多的不超過 \(n\) 的數是多少？

解題思路

思路同上，只不過要改改 dfs 的返回條件。注意這樣的題目的數據範圍，32 位整數可能溢出。

參考代碼

#include <cstdio>
#include <iostream>

int p[16] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53};
unsigned long long n;
unsigned long long ans,
    ans_num;  // ans 为 n 以内的最大反素数（会持续更新），ans_sum 为
              // ans的因子数。

// depth: 当前在枚举第几个素数
// temp: 当前因子数量为 num的时候的数值
// num: 当前因子数
// up：上一个素数的幂，这次应该小于等于这个幂次嘛
void dfs(int depth, unsigned long long temp, unsigned long long num, int up) {
  if (depth >= 16 || temp > n) return;
  if (num > ans_num) {  // 更新答案
    ans = temp;
    ans_num = num;
  }
  if (num == ans_num && ans > temp) ans = temp;  // 更新答案
  for (int i = 1; i <= up; i++) {
    if (temp * p[depth] > n)
      break;  // 剪枝：如果加一个这个乘数的结果比ans要大，则必不是最佳方案
    dfs(depth + 1, temp *= p[depth], num * (i + 1),
        i);  // 取一个该乘数，进行对下一个乘数的搜索
  }
  return;
}

int main() {
  while (scanf("%llu", &n) != EOF) {
    ans_num = 0;
    dfs(0, 1, 1, 60);
    printf("%llu\n", ans);
  }
  return 0;
}

參考資料與註釋

Rui-Juan Jing, Marc Moreno-Maza, Delaram Talaashrafi, "Complexity Estimates for Fourier-Motzkin Elimination", Journal of Functional Programming 16:2 (2006) pp 197-217.
數論部分第一節：素數與素性測試
Miller-Rabin 與 Pollard-Rho 學習筆記 - Bill Yang's Blog
Primality test - Wikipedia
桃子的算法筆記——反素數詳解（acm/OI）
The Rabin-Miller Primality Test
Bach, Eric , "Explicit bounds for primality testing and related problems", Mathematics of Computation, 55:191 (1990) pp 355–380.
Deterministic variant of the Miller-Rabin primality test
Fermat pseudoprime - Wikipedia
Carmichael number - Wikipedia
Carmichael function - Wikipedia
Carmichael Number -- from Wolfram MathWorld
Carmichael's Lambda Function | Brilliant Math & Science Wiki

Korselt, A. R. (1899). "Problème chinois".L'Intermédiaire des Mathématiciens.6: 142–143. ↩
W. R. Alford; Andrew Granville; Carl Pomerance (1994). "There are Infinitely Many Carmichael Numbers".Annals of Mathematics. 140 (3): 703–722. ↩
Erdős, P. (1956). "On pseudoprimes and Carmichael numbers".Publ. Math. Debrecen. 4 (3–4): 201–206. ↩
PINCH, Richard GE. The Carmichael numbers up to \({10}^{20}\). ↩

本页面最近更新：，更新历史
发现错误？想一起完善？在 GitHub 上编辑此页！
本页面贡献者：OI-wiki
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用