梅森素数

梅森素数 (Mersenne prime) 是人类目前已知的最大的一类素数. 它是指形如 `M_p = 2^p - 1` 的素数. 若 `p = a b`, 则 `2^(a b) - 1` 可以因式分解. 因此 `2^p - 1` 是素数的前提条件是 `p` 为一素数. 前几个梅森素数的 `p` 值是 2, 3, 5, 7, 13, 17, 19, 31, 61, 89, 107, 127... 晚年的 Euler 在双目失明的情况下通过心算验证 `M_31 = 2147483647` 为素数.

完全数 (perfect number) 是指满足方程 `sigma(n) = 2n` 的正整数. 其中 `sigma(n)` 表示 `n` 的所有正因子之和. 前几个完全数是 6, 28, 496, 8128, 33550336... 它们正好是偶数. 奇完全数是否存在仍是一个未解之谜.

梅森素数与偶完全数的一一对应 若 `M_p` 为梅森素数, 则相应的三角形数 `M_p(M_p+1)//2` 是完全数. 反之, 偶完全数必然具有这种形式.

  1. [古希腊, Euclid] 考虑 `2^(p-1) * (2^p-1)` 的全体正因子 `1, 2, cdots, 2^(p-1)`,
    `M_p, 2 M_p, cdots, 2^(p-1) M_p`.
    求和即得结论.
  2. [Leonhard Euler] 下证偶完全数必然形如 `M_p(M_p+1)//2`. 设这个偶完全数是 `2^a b`, 其中 `a` 是正整数, `b` 是奇数. 由于 `sigma` 是积性函数, 有 `2^(a+1) b` `= sigma(2^a b)` `= sigma(2^a) sigma(b)` `= (2^(a+1) - 1) sigma(b)`. 记 `n = a+1`, 则 `2^n b = (2^n - 1) sigma(b)`. 这指出 `2^n - 1` 整除 `b`, 且 `2^n` 整除 `sigma(b)`. 设 `sigma(b) = 2^n c` 则 `b = (2^n - 1) c`. 假如 `c gt 1`, 显然 `2^n - 1` 也大于 `1`, 则 `b` 至少有三个不同因子 `1, c, b`, 于是 `sigma(b) ge 1 + c + b` `= 1 + c + (2^n - 1) c` `= 1 + 2^n c` `= 1 + sigma(b)`, 矛盾. 因此 `c = 1`, 即 `b = 2^n - 1`, `sigma(b) = 2^n`. 这又推出 `b` 是素数.

设 `p` 为奇素数, 则 `2^p - 1` 的所有因子都形如 `2n p+1`, `n` 是正整数.

由于两个形如 `2n p + 1` 的数相乘, 结果仍为同一形式, 我们只需对素数作出证明. 令 `q` 是 `2^p - 1` 的素因子, 则 `q` 是奇素数. 由费马小定理, `q | 2^(q-1) - 1`. 利用公式 `gcd(2^a - 1, 2^b - 1) = 2^(gcd(a,b)) - 1`, 我们有 `2^(gcd(p, q-1)) - 1` `ge q gt 1`, 即 `gcd(p, q-1) gt 1`, 即 `p | q-1`. 又 `q-1` 是偶数, 所以存在正整数 `n` 使得 `q = 2n p + 1`.

Lucas-Lehmer 算法 用于验证 `M_p` 是否为梅森素数. 设 `p` 是奇素数, `M_p = 2^p - 1`. 令 `r_1 = 4`, `r_n = (r_(n-1)^2 - 2) mod M_p`, (取最小非负剩余) 则 `M_p` 是素数当且仅当 `r_(p-1) = 0`.

# 该程序在个人电脑上验证第 27 个梅森素数 2**44497-1 用时两分钟
# 截至 2026 年, 人类已发现 52 个梅森素数
def test_mersenne(p):
    if p == 2:
        return True
    mp = 2**p-1
    r = 4
    for i in range(p-2):
        r = pow(r, 2, mp)
        r = (r-2) % mp
    return r == 0

素数定理

关于素数分布的一个重要结论是素数定理: `pi(x) ~ "li"(x) ~ x/(ln x)`, `quad x to oo`. 其中 `pi(x)` 是 `x` 以内 (即 `le x`) 的素数个数, 而 `"li"(x) = int_0^x dt/(ln t)` `= lim_(epsi to 0^+) int_0^(1-epsi) + int_(1+epsi)^x`. 这意味着在 `x` 附近的素数密度约为 `1//ln x`. 素数定理的证明比较复杂, 我们只证较弱的结论 `pi(x) = Theta(x/(ln x))`. 换言之, `x` 充分大时, 存在正的常数 `c_1, c_2` 使 `c_1 x/(ln x) le pi(x) le c_2 x/(ln x)`. 此即著名的 Chebyshev 不等式.

由素数定理推出 `p_n ~ n ln n`, `p_n` 是第 `n` 个素数.

[群友@zz] 由 `pi(x) ~ x/(ln x)` 取对数得 `ln pi(x) ~ ln x - ln ln x`. 于是利用 `n = pi(p_n)` 有 `p_n/(n ln n)` `= p_n/(pi(p_n)) * 1/ln(pi(p_n))` `~ ln p_n * 1/(ln p_n - ln ln p_n)` `~ 1`.

素数之比在 `RR^+` 中稠密.

[知乎@用户pKD90B] 任取 `alpha gt 0`, 令正数列 `alpha_j to alpha` `(j ge 2)`, 并定义 `m_j := |__j // ln j__|`, `quad n_j := |__j alpha_j//ln j__|`, 下证 `lim_(j to oo) p_(n_j)//p_(m_j) = alpha`. 由于 `p_n ~ n ln n`, 即证 `(n_j ln n_j)/(m_j ln m_j) to alpha`. 事实上, 由 `(j/(ln j) - 1) ln(j/(ln j) - 1)` `lt m_j ln m_j` `le j/(ln j)(ln j - ln ln j)` 知道, `m_j ln m_j ~ j`. 类似地, 由 `((j alpha_j)/(ln j) - 1) ln ((j alpha_j)/(ln j) - 1)` `lt n_j ln n_j` `le (j alpha_j)/(ln j) (ln j + ln alpha_j - ln ln j)` 知道, `n_j ln n_j ~ j alpha_j`.

筛法

用 Eratosthenes 筛法求 100 以内的素数个数.

考虑 `sqrt 100 = 10` 以内的全体素数 2, 3, 5, 7. 对于 `10 lt n le 100`, 只要 `n` 与 2, 3, 5, 7 互素, 就能保证为素数. 记 `A = {1, cdots, 100}`, `P = 2 * 3 * 5 * 7`, 考查集合 `A` 中与 `P` 互素的整数个数, 称为筛函数: `S(A, P) = \#{a in A: (a, P) = 1}` `= sum_(a in A, (a,P) = 1) 1`. 在我们的例子中, 100 以内与 2, 3, 5, 7 互素的整数包括 10 到 100 的所有素数, 以及整数 1. 于是 `S(A, P) = pi(100) - pi(10) + 1`. 接下来利用容斥原理计算筛函数. 记 `A_d` 为 `A` 中所有 `d` 的倍数, 则 `|A_d| = |__100//d__|`. `S(A, P)` `= |A| - |A_2| - |A_3| - |A_5| - |A_7|`
`+ |A_(2*3)| + |A_(2*5)| + |A_(2*7)| + |A_(3*5)| + |A_(3*7)| + |A_(5*7)|`
`- |A_(2*3*5)| - |A_(2*3*7)| - |A_(2*5*7)| - |A_(3*5*7)|` `+ |A_(2*3*5*7)|`
`= 22`.
因此 `pi(100) = S(A, P) + pi(10) - 1 = 22 + 4 - 1 = 25`.

使用容斥原理求筛函数 设 `A` 为有限个整数的集合, `A_d` 为 `A` 中所有 `d` 的倍数. 又正整数 `P` 的所有素因子为 `p_1, cdots, p_s`, 那么根据容斥原理, `S(A, P) := sum_(a in A, (a,P) = 1) 1` `= sum_(r=0)^s (-1)^r sum_(i_1 lt cdots lt i_r) |A_(p_(i_1) cdots p_(i_r))|`. 使用上一章的 Möbius 函数, 结论简化为 `S(A, P) = sum_(d | P) mu(d) |A_d|`. 特别当 `A = {1, cdots, n}`, `P = n` 时, 得到 Euler 函数的表达式 `varphi(n) = sum_(d | n) mu(d) n//d` `= n prod_(p | n) sum_(p^a | n) mu(p^a)//p^a` `= n prod_(p | n) (1 - 1//p)`.

利用公式 `sum_(d | P) mu(d) = { 1, if P = 1; 0, if P gt 1 :}` `sum_(a in A, (a,P) = 1) 1` `= sum_(a in A) sum_(d | (a,P)) mu(d)` `= sum_(d | P) mu(d) sum_(a in A, d | a) 1` `= sum_(d | P) mu(d) |A_d|`.

粗略估算, 以 `le a` 的素数去筛 `n` 以内的正整数时, 每个素数 `p` 大约筛掉 `1//p` 的选项, 因此通过筛选的数大约有这么多: `n prod_(p le a) (1-1//p)`.

    筛法求 `pi(n)` 用前 `a` 个素数去筛 `x` 以内的整数, 考虑筛函数 `varphi(x, a) := \#{ n le x: gcd(n, p_1 cdots p_a) = 1}`.
  1. 递推公式 `varphi(x, a) = varphi(x, a-1) - varphi(x//p_a, a-1)`.
  2. Legendre 等式 `pi(x) = varphi(x, a) + a - 1`, 其中 `a = pi(sqrt x)`.
  3. Lehmer 公式 `pi(x) = varphi(x, a) + a - 1 - P_2`, 其中 `a = pi(root 3 x)`. `P_2` 是区间 `(a, x]` 中通过筛选的合数的数量, 这些合数必定拥有两个素因子, 且最小素因子大于 `a`. 事实上 `P_2 = sum_(root 3 x lt p le sqrt x) (pi(x//p) - pi(p) + 1)`, `quad p` 为素数.
  4. 实测 Legendre 公式在 `10^8` 量级就会超过 python 的递归深度,但 Lehmer 不会.
  1. `x` 以内的整数要与前 `a` 个素数都互素, 首先要与前 `a-1` 个素数互素. 然后扣除掉漏网之鱼, 即那些与前 `a-1` 个素数互素, 但满足 `p_a | n` 的整数 `n`. 这当且仅当 `n//p_a` 与前 `a-1` 个素数互素.
  2. 这是因为 `x` 以内的合数最小素因子不超过 `sqrt x`. 换言之只要一个大于 `sqrt x` 的数通过了筛选, 它必为素数.
  3. `x` 以内的合数, 若它的素因子个数 `ge 3` (按重数计算), 则它的最小素因子不超过 `root 3 x`. 换言之一个大于 `root 3 x` 的数如果通过筛选, 它要么是素数, 要么是一个拥有两个因子的合数, 且最小素因子大于 `root 3 x`, 小于等于 `sqrt x`. 记这个合数为 `p q`, 两个素因子满足 `p le q le x//p`. 于是因子 `q` 的可能取法有 `pi(x//p) - pi(p) + 1` 种, 求和就得到 `P_2` 表达式.

Chebyshev 不等式

    Bertrand-Chebyshev 定理 [知乎@Zephyr, 数学天书中的证明]
  1. 对任意正整数 `n` 有 `pi(n) lt pi(2n)`, 换言之存在素数 `p` 满足 `n lt p le 2n`. 当 `n ge 2` 时, `2n` 不可能是素数, 不等式严格成立.
  2. 由 1. 立即得到较松的估计: `pi(2^n) ge n`, `p_n le 2^n`.
  1. 先证 `n` 以内的素数乘积上界: `prod_(p le n) p le 4^(n-1)`. `n le 2` 时直接验证成立. `n ge 3` 时, 只需考虑奇素数的情形. 设 `n = 2m+1`, 若结论对一切整数 `2 le n le 2m` 都成立, 则 `prod_(p le 2m+1) p` `= prod_(p le m+1) p prod_(m+1 lt p le 2m+1) p` `le 4^m (2m+1;m)` `= 4^m ((2m;m) + (2m;m-1))` `le 4^m (1+1)^(2m)` `= 4^(2m)`. 故结论成立. 注意把 `n` 换成 `ge 1` 的实数, 不等式仍成立.
  2. 估算每个素因子 `p` 在 `N = (2n;n)` 中的次数. 事实上, 取 `r ge 0` 使得 `p^r le 2n lt p^(r+1)`, 则 `v_p(N)` `= v_p((2n)!) - 2 v_p(n!)` `= sum_(k=1)^r (floor((2n)/p^k) - 2floor(n/p^k))` `le r`. 这推出 `p^(v_p(N)) le p^r le 2n`.
    取 `r = 1` 得到: `p gt sqrt(2n)` 时, `v_p(N) le 1`.
    又设 `n ge 5`, 此时有 `2 lt sqrt(2n) lt 2n//3`, 对于 `2n//3 lt p le n`, 直接观察 `N = (2n)!/(n! n!)` 的分子分母可知 `v_p(N) = 0`.
  3. 现在考虑 `N` 的素因子分解. 当 `n ge 5` 时, `N` 有上界 `N = prod_(p le 2n) p^(v_p(N))` `le prod_(2 le p le sqrt(2n)) 2n prod_(sqrt(2n) lt p le 2n//3) p prod_(n lt p le 2n) p`. 另一方面, 不等式 `2n (2n;n) ge (1+1)^(2n)` 给出 `N` 的下界: `N ge 4^n//2n`. 假设不存在素数 `p` 满足 `n lt p le 2n`, 由以上两式, 结合 1. 的结论推出 `4^n le 2n * N` `le (2n)^(sqrt(2n)) prod_(sqrt(2n) lt p le 2n//3) p` `le (2n)^(sqrt(2n)) 4^(2n//3)`. 于是 `4^(n//3) le (2n)^(sqrt(2n))`, 即 `sqrt(2n) le 3 op(log_2)(2n)`, 这推出 `n le 426`. 但我们可以直接验证 Bertrand 结论对 `n lt 500` 成立, 例如, 下面相邻两项素数的比值小于 2: 2, 3, 5, 7, 13, 23, 43, 83, 163, 317, 631. 因此假设不正确, 证毕.

Chebyshev 不等式 记 `pi(x)` 是不超过 `x` 的素数个数, `p_n` 为第 `n` 个素数. 我们有 `pi(2m) ge m // log_2(2m)`, `quad AA m ge 1`,
`pi(2m) - pi(m) lt 2m // log_2 m`, `quad AA m ge 1`,
`1/3 * x/(log_2 x)` `lt pi(x)` `lt 6 * x/(log_2 x)`, `quad AA x ge 2`,
`1/6 * n log_2 n` `lt p_n` `lt 8 * n log_2 n`, `quad AA n ge 2`.

    先证前两个不等式. 我们从 `M = (2 m; m)` 的估计入手, `m` 为正整数.
  1. 联系第一章 Kummer 定理的相关内容, 记 `v_p(n)` 为素数 `p` 在 `n` 中的次数, 有 `v_p(n!) = sum_(i ge 1) |__n/p^i__|`. 取对数有 `ln M` `= ln((2m)!) - 2 ln(m!)`
    `= sum_(p lt 2m) v_p((2m)!) ln p` `- 2 sum_(p le m) v_p(m!) ln p`
    `= sum_(p le m) [v_p((2m)!) - 2 v_p(m!)] ln p` `+ sum_(m lt p lt 2m) v_p((2m)!) ln p`
    `= A + B`.
  2. 对 `y gt 0` 有 `0 le |__2y__| - 2|__y__| le 1`, 于是 `0 le (|__(2m)/p^i__| - 2|__m/p^i__|) le 1`, `quad i ge 1`. 上式求和得到 `A` 的系数的估计 `0 le v_p((2m)!) - 2 v_p(m!)` `le sum_(p^i | 2m) 1` `= |__ln(2m)/ln p__|`. 另一方面, `m lt p lt 2m` 时, 显然 `p` 在 `(2m)!` 中的次数为 `1`, 即 `B` 的系数 `v_p((2m)!) = 1`. 于是 `sum_(m lt p lt 2m) ln p` `le ln M` `le sum_(p lt 2m) |__ln(2m)/ln p__| ln p`. 因而 `(pi(2m)-pi(m)) ln m` `le ln M` `le pi(2m) ln(2m)`.
  3. 通过直接估计 `M` 的上下界得到 `M = (2m)/m * (2m-1)/(m-1) cdots (m+1)/1 ge 2^m`;
    `M = (2m; m) lt (1+1)^(2m) = 2^(2m)`.
    于是 `pi(2m) ln(2m) ge m ln 2`,
    `(pi(2m)-pi(m))ln m lt 2m ln 2`.
    应用前两个不等式证 `pi(x)` 的上下界.
  1. 当 `x ge 6` 时, 取 `m = |__x//2__| gt 2`, 此时成立 `2m lt x lt 3m`, 因此 `pi(x) log_2 x ge pi(2m) log_2(2m)` `ge m` `gt x/3`. 直接验算知上式对 `2 le x lt 6` 也成立; 这证明了不等式的左半部分.
  2. 当 `m = 2^k` 时, `(pi(2^(k+1)) - pi(2^k)) k lt 2^(k+1)`, 又显然 `pi(2^(k+1)) le 2^k`, 两式相加得 `(k+1)pi(2^(k+1)) - k pi(2^k) lt 3 * 2^k`. 令 `k` 从 0 到 `h-1` 求和有 `h pi(2^h) lt 3 * 2^h`. 设 `2^(h-1) lt x le 2^h`, 则 `pi(x) le pi(2^h)` `lt 3 * 2^h//h` `lt 3 * (2 x)/(log_2 x)`. 这证明了不等式的右半部分.
    再来证 `p_n` 的上下界.
  1. 在 `pi(x)` 的估计中代入 `x = p_n`. 注意 `pi(p_n) = n` 和 `p_n gt n`, 有 `n lt 6 p_n/(log_2 p_n)` `lt 6 p_n/(log_2 n)`, 这证明了左半不等式.
  2. `n = 2` 时直接验证右半不等式成立. `n ge 3` 时 `p_n` 为奇素数. 在不等式 `pi(2m) ge m // log_2 (2m)` 中取 `2m = p_n + 1`, 得 `n ln(p_n+1) ge (p_n+1)/2 * ln 2`. 整理得 `ln(p_n+1) le ln(2n//ln 2) + ln ln(p_n+1)`. 注意 `n ge 3 gt 2//ln 2`, 并利用不等式 `ln x lt x//2` (`x gt 0`), 有 `ln(p_n+1) le 2 ln n + 1/2 ln(p_n+1)`. 上式结合 `(p_n+1)/(2n // ln 2)` `le ln(p_n+1)` `lt 4 ln n`. 这证明了右半不等式.
    Chebyshev 函数 两个 Chebyshev 函数定义为
  1. `vartheta(x) = sum_(p le x) ln p`, 其中 `p` 为素数;
  2. `psi(x) = sum_(p^a le x) ln p`, 其中 `p` 为素数, `a` 为正整数.
  3. 若定义 Mangoldt 函数为 `Lambda(n) = { ln p, if n = p^a; 0, "otherwise" :}`, 则 `psi(x) = sum_(n le x) Lambda(n)`.

Mangoldt 函数的和函数是对数函数: `sum_(d | n) Lambda(d) = ln n`.

设 `p` 在 `n` 中的次数为 `a_p`, 利用算术基本定理, 左边等于 `sum_(p^a | n) ln p` `= a_p sum_(p | n) ln p` `= ln(prod_(p | n) p^(a_p))` `= ln n`.

`x ge 1` 时, `sum_(n le x) psi(x // n) = ln(|__x__|!)`.

由定义, 左边等于 `sum_(n le x) sum_(k le x//n) Lambda(k)` `= sum_(n le x) sum_(n k le x) Lambda(k)` `==^(n k = m) sum_(m le x) sum_(n | m) Lambda(m/n)` `= sum_(m le x) ln m` 等于右边.

    Chebyshev 函数是 `ln x pi(x)` 的良好替代, 事实上存在常数 `c_gamma gt 0` 使得
  1. `(ln x - c_gamma) pi(x) lt vartheta(x) lt ln x pi(x)`;
  2. `vartheta(x) le psi(x) le vartheta(x) + sqrt x ln x`.
    以下三个命题等价:
  1. 存在常数 `d_1, d_2 gt 0` 使得 `d_1 x // ln x lt pi(x) lt d_2 x // ln x`;
  2. 存在常数 `d_3, d_4 gt 0` 使得 `d_3 x lt vartheta(x) lt d_4 x`;
  3. 存在常数 `d_5, d_6 gt 0` 使得 `d_5 x lt psi(x) lt d_6 x`;
  4. 以下三个命题等价: `x to oo` 时
  5. `pi(x) ~ x // ln x`;
  6. `vartheta(x) ~ x`;
  7. `psi(x) ~ x`.