P3383 【模板】线性筛素数 题解

· · 题解

这是一篇 Python 题解。

{\color{#00AACD}\textbf{埃拉托斯特尼筛法（埃氏筛）}}

{\color{#00CD00}\text{算法介绍}}

根据质数的定义，任何一个质数的 x 倍（x > 1）都是合数。于是我们可以从小到大考虑每一个数，并将这个数的所有大于自身倍数都标记为合数。运行结束后，所有未被标记的数就是质数了。

算法的流程如下：

  1. 创建一个初始时全为 1 的标记数组 isprime，并将 0 和 1 标记为非质数。
  2. 从小到大遍历 2\sim n 的每个整数：
    • 如果这个数已经被标记为合数，则将其跳过；
    • 否则，遍历这个数的所有大于自身的倍数，将它们都标记为合数。
  3. 最后，所有未被标记的数即为质数。

以上即为埃拉托斯特尼筛法（埃氏筛）。埃氏筛的时间复杂度为 O(n\log\log n)，证明较复杂，可参见 OI - Wiki。

虽然埃氏筛的时间复杂度是 O(n\log\log n)，但是由于其常数较小，且有多种方式进行优化，很多时候跑得比 O(n) 的欧拉筛还要快。

{\color{#00CD00}\text{算法优化}}

  1. 要找出 1\sim n 的所有质数，只需要对 1\sim\lfloor\sqrt n\rfloor 的质数进行筛选就够了。正确性显然。
  2. 除了 2 以外的所有偶数都是合数，因此只需要筛奇数就行了，运算量直接减少了一半。
  3. 对于一个质数 p，只需从它的 p 倍开始标记即可，因为 2\sim p-1 倍已经被前面的质数标记了。
  4. 标记数组 isprime 可以用 bool 类型存储，可以减少内存占用。在 Python 中，可以用 NumPy 库创建数组。

{\color{#00CD00}\text{代码实现}}

import numpy as np # 使用 NumPy 库提供高效的列表操作
from sys import stdin, stdout # 使用 sys.stdin/stdout 代替缓慢的 input/print

n, q = map(int, stdin.readline().split(' '))
isprime = np.ones(n + 1, dtype = bool) # 创建初始全为 1，类型为 bool 的列表
isprime[0] = isprime[1] = isprime[4::2] = False # 将 0、1 和大于 2 的偶数标记为非质数

for i in range(3, int(n**0.5) + 1, 2): # 只筛 sqrt(n) 以内的奇数
    if isprime[i]: isprime[i*i::i+i] = False # 从 i*i 开始筛，步长为 2i

prime = np.nonzero(isprime)[0] # 使用 np.nonzero() 找出 isprime 中所有非 0 的索引

for k in stdin: stdout.write(str(prime[int(k) - 1]) + '\n')

{\color{#00AACD}\textbf{欧拉筛法（线性筛）}}

{\color{#00CD00}\text{算法介绍}}

注意到在埃氏筛中，一个数可能会被多个质数标记（如 12 = 2\times 6 = 3\times 4，既可以被 2 标记也可以被 3 标记）。考虑如何让一个合数只被其最小的质因数标记。设计如下算法：

  1. 创建一个初始时全为 1 的标记数组 isprime 和一个初始时为的质数列表 prime。
  2. 从小到大遍历 2\sim n 的每个数，设当前的数为 i，如果 i 未被标记，就将其加入质数列表中。
    • 对于每个 i，遍历当前质数列表中的每一个质数，设当前的质数为 p，将 i\times p 标记为合数。
    • 如果 i\bmod p=0，说明 i 之前已经被 p 标记过了，所以 i 乘上其他质数的结果也一定会被 p 标记，因此直接退出即可。
  3. 最后，得到的质数列表 prime 即为 1\sim n 的所有质数。

以上即为欧拉筛法（线性筛）。可以发现，在线性筛中，每个合数都只会被其最小的质因数标记 1 次，因此时间复杂度为 O(n)。

{\color{#00CD00}\text{代码实现}}

在本题中，线性筛做法需要使用 PyPy3 提交才能勉强通过。洛谷的 PyPy3 不提供 NumPy 库支持，故使用 bytearray 代替。

from sys import stdin, stdout # 使用 sys.stdin/stdout 代替缓慢的 input/print

n, q = map(int, stdin.readline().split(' '))
isprime = bytearray([1]) * (n + 1) # 创建长度为 n+1, 初始全为 1 的 bytearray
prime = [] # 创建初始为空的质数列表 prime

for i in range(2, n + 1): # 遍历 2~n 的每一个数
    if isprime[i]: prime.append(i) # 如果当前的数未被标记，将其加入质数列表
    for p in prime: # 遍历当前质数列表中的每一个质数
        if i * p > n: break # 若超出范围直接退出
        isprime[i * p] = 0 # 将 i*p 标记为合数
        if i % p == 0: break # 若 i%p==0 直接退出

for k in stdin: stdout.write(str(prime[int(k) - 1]) + '\n')