列出N以下所有素数的最快方法

Question

列出N以下所有素数的最快方法这是我能提出的最佳算法。def get_primes(n):    numbers = set(range(n, 1, -1))    primes = []    while numbers:        p = numbers.pop()        primes.append(p)        numbers.difference_update(set(range(p*2, n+1, p)))    return primes>>> timeit.Timer(stmt='get_primes.get_primes(1000000)', setup='import   get_primes').timeit(1)1.1499958793645562可以做得更快吗？此代码有一个缺陷：由于numbers是无序集，因此无法保证numbers.pop()从集中删除最小数字。然而，它对某些输入数字起作用（至少对我而言）：>>> sum(get_primes(2000000))142913828922L#That's the correct sum of all numbers below 2 million>>> 529 in get_primes(1000)False>>> 529 in get_primes(530)True

Answer 1

更快，更内存的纯Python代码：

def primes(n):
    """ Returns  a list of primes < n """
    sieve = [True] * n    for i in range(3,int(n**0.5)+1,2):
        if sieve[i]:
            sieve[i*i::2*i]=[False]*((n-i*i-1)//(2*i)+1)
    return [2] + [i for i in range(3,n,2) if sieve[i]]

或者从半筛开始

def primes1(n):
    """ Returns  a list of primes < n """
    sieve = [True] * (n//2)
    for i in range(3,int(n**0.5)+1,2):
        if sieve[i//2]:
            sieve[i*i//2::i] = [False] * ((n-i*i-1)//(2*i)+1)
    return [2] + [2*i+1 for i in range(1,n//2) if sieve[i]]

更快，更内存的numpy代码：

import numpydef primesfrom3to(n):
    """ Returns a array of primes, 3 <= p < n """
    sieve = numpy.ones(n//2, dtype=numpy.bool)
    for i in range(3,int(n**0.5)+1,2):
        if sieve[i//2]:
            sieve[i*i//2::i] = False
    return 2*numpy.nonzero(sieve)[0][1::]+1

从筛子的三分之一开始变化更快：

import numpydef primesfrom2to(n):
    """ Input n>=6, Returns a array of primes, 2 <= p < n """
    sieve = numpy.ones(n//3 + (n%6==2), dtype=numpy.bool)
    for i in range(1,int(n**0.5)//3+1):
        if sieve[i]:
            k=3*i+1|1
            sieve[       k*k//3     ::2*k] = False
            sieve[k*(k-2*(i&1)+4)//3::2*k] = False
    return numpy.r_[2,3,((3*numpy.nonzero(sieve)[0][1:]+1)|1)]

上述代码的（难以编码）纯python版本将是：

def primes2(n):
    """ Input n>=6, Returns a list of primes, 2 <= p < n """
    n, correction = n-n%6+6, 2-(n%6>1)
    sieve = [True] * (n//3)
    for i in range(1,int(n**0.5)//3+1):
      if sieve[i]:
        k=3*i+1|1
        sieve[      k*k//3      ::2*k] = [False] * ((n//6-k*k//6-1)//k+1)
        sieve[k*(k-2*(i&1)+4)//3::2*k] = [False] * ((n//6-k*(k-2*(i&1)+4)//6-1)//k+1)
    return [2,3] + [3*i+1|1 for i in range(1,n//3-correction) if sieve[i]]

不幸的是，pure-python不采用更简单，更快速的numpy方式进行赋值，并且len()在循环内调用[False]*len(sieve[((k*k)//3)::2*k])太慢。所以我不得不即兴改正输入（避免更多的数学运算）并做一些极端（和痛苦）的数学魔术。

就个人而言，我认为numpy（这是如此广泛使用）不是Python标准库的一部分是一种耻辱，并且Python开发人员似乎完全忽略了语法和速度的改进。