如何高效寻找素数

素数的定义看起来很简单，如果一个数如果只能被 1 和它本身整除，那么这个数就是素数。

不要觉得素数的定义简单，恐怕没多少人真的能把素数相关的算法写得高效。比如让你写这样一个函数：

  1. // 返回区间 [2, n) 中有几个素数 
  2. int countPrimes(int n)
  4. // 比如 countPrimes(10) 返回 4
  5. // 因为 2,3,5,7 是素数

你会如何写这个函数？我想大家应该会这样写：

  1. int countPrimes(int n) {
  2.     int count = 0;
  3.     for (int i = 2; i < n; i++)
  4.         if (isPrim(i)) count++;
  5.     return count;
  6. }
  8. // 判断整数 n 是否是素数
  9. boolean isPrime(int n) {
  10.     for (int i = 2; i < n; i++)
  11.         if (n % i == 0)
  12.             // 有其他整除因子
  13.             return false;
  14.     return true;
  15. }

这样写的话时间复杂度 O(n^2)，问题很大。首先你用 isPrime 函数来辅助的思路就不够高效；而且就算你要用 isPrime 函数，这样写算法也是存在计算冗余的。

先来简单说下如果你要判断一个数是不是素数，应该如何写算法。只需稍微修改一下上面的 isPrim 代码中的 for 循环条件：

  1. boolean isPrime(int n) {
  2.     for (int i = 2; i * i <= n; i++)
  3.         ...
  4. }

换句话说，i 不需要遍历到 n，而只需要到 sqrt(n) 即可。为什么呢，我们举个例子，假设 n = 12。

  1. 12 = 2 × 6
  2. 12 = 3 × 4
  3. 12 = sqrt(12) × sqrt(12)
  4. 12 = 4 × 3
  5. 12 = 6 × 2

可以看到，后两个乘积就是前面两个反过来，反转临界点就在 sqrt(n)。

换句话说，如果在 [2,sqrt(n)] 这个区间之内没有发现可整除因子，就可以直接断定 n 是素数了，因为在区间 [sqrt(n),n] 也一定不会发现可整除因子。

现在，isPrime 函数的时间复杂度降为 O(sqrt(N))，但是我们实现 countPrimes 函数其实并不需要这个函数，以上只是希望读者明白 sqrt(n) 的含义，因为等会还会用到。

高效实现 countPrimes

高效解决这个问题的核心思路是和上面的常规思路反着来：

首先从 2 开始，我们知道 2 是一个素数，那么 2 × 2 = 4, 3 × 2 = 6, 4 × 2 = 8… 都不可能是素数了。

然后我们发现 3 也是素数，那么 3 × 2 = 6, 3 × 3 = 9, 3 × 4 = 12… 也都不可能是素数了。

看到这里，你是否有点明白这个排除法的逻辑了呢？先看我们的第一版代码：

  1. int countPrimes(int n) {
  2.     boolean[] isPrim = new boolean[n];
  3.     // 将数组都初始化为 true
  4.     Arrays.fill(isPrim, true);
  6.     for (int i = 2; i < n; i++) 
  7.         if (isPrim[i]) 
  8.             // i 的倍数不可能是素数了
  9.             for (int j = 2 * i; j < n; j += i) 
  10.                     isPrim[j] = false;
  12.     int count = 0;
  13.     for (int i = 2; i < n; i++)
  14.         if (isPrim[i]) count++;
  16.     return count;
  17. }

如果上面这段代码你能够理解，那么你已经掌握了整体思路，但是还有两个细微的地方可以优化。

首先，回想刚才判断一个数是否是素数的 isPrime 函数，由于因子的对称性，其中的 for 循环只需要遍历 [2,sqrt(n)] 就够了。这里也是类似的，我们外层的 for 循环也只需要遍历到 sqrt(n)：

  1. for (int i = 2; i * i < n; i++) 
  2.     if (isPrim[i]) 
  3.         ...

除此之外，很难注意到内层的 for 循环也可以优化。我们之前的做法是：

  1. for (int j = 2 * i; j < n; j += i) 
  2.     isPrim[j] = false;

这样可以把 i 的整数倍都标记为 false，但是仍然存在计算冗余。

比如 n = 25，i = 4 时算法会标记 4 × 2 = 8，4 × 3 = 12 等等数字，但是这两个数字已经被 i = 2 和 i = 3 的 2 × 4 和 3 × 4 标记了。

我们可以稍微优化一下，让 j 从 i 的平方开始遍历，而不是从 2 * i 开始：

  1. for (int j = i * i; j < n; j += i) 
  2.     isPrim[j] = false;

这样，素数计数的算法就高效实现了，其实这个算法有一个名字，叫做 Sieve of Eratosthenes。看下完整的最终代码：

  1. int countPrimes(int n) {
  2.     boolean[] isPrim = new boolean[n];
  3.     Arrays.fill(isPrim, true);
  4.     for (int i = 2; i * i < n; i++) 
  5.         if (isPrim[i]) 
  6.             for (int j = i * i; j < n; j += i) 
  7.                 isPrim[j] = false;
  9.     int count = 0;
  10.     for (int i = 2; i < n; i++)
  11.         if (isPrim[i]) count++;
  13.     return count;
  14. }

该算法的时间复杂度比较难算，显然时间跟这两个嵌套的 for 循环有关，其操作数应该是：

n/2 + n/3 + n/5 + n/7 + … = n × (1/2 + 1/3 + 1/5 + 1/7…)

括号中是素数的倒数。其最终结果是 O(N * loglogN)，有兴趣的读者可以查一下该算法的时间复杂度证明。

以上就是素数算法相关的全部内容。怎么样，是不是看似简单的问题却有不少细节可以打磨呀？

致力于把算法讲清楚！欢迎关注我的微信公众号 labuladong，查看更多通俗易懂的文章：

labuladong

ruicore 提供 Python3 代码：

  1. class Solution:
  2.     def countPrimes(self, n: int):
  3.         # 前 2 个数不是素数
  4.         if n < 3:
  5.             return 0
  6.         # isprime 数组
  7.         isprime = [1] * n
  9.         for i in range(2, int(n ** 0.5) + 1):
  10.             if isprime[i]:
  11.                 # 从 i*i（包括自己） 到 n（不包括n），每步增加 i 的所有数的个数
  12.                 tmp = ((n - 1 - i * i) // i + 1)
  13.                 # 这种方式赋值比用 for 循环更高效
  14.                 isprime[i * i: n: i] = [0] * tmp
  16.         # 前 2 个数不是素数，去掉
  17.         count = sum(isprime[2:])
  18.         return count