这是我能想到的最好的算法。

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

当前回答

我知道比赛已经结束好几年了。...

尽管如此,这是我对纯python质数筛子的建议,基于在向前处理筛子时使用适当的步骤省略2、3和5的倍数。尽管如此,在N<10^9时,它实际上比@Robert William Hanks的优解rwh_primes2和rwh_primes1要慢。通过使用大于1.5* 10^8的ctypes.c_ushort筛分数组,可以在某种程度上适应内存限制。

10^6

$ python -mtimeit -s"import primeSieveSpeedComp" "primeSieveSpeedComp. primesieveseq (1000000)" 10个循环,最好的3:46.7毫秒每循环

import primeSieveSpeedComp (primeSieveSpeedComp) “primeSieveSpeedComp.rwh_primes1(1000000)”10个循环,最好的3:43.2 每回路Msec $ python -m timeit -s"import primeSieveSpeedComp" “primeSieveSpeedComp.rwh_primes2(1000000)”10圈,最好成绩是3:34.5 每回路Msec

10^7

$ python -mtimeit -s"import primeSieveSpeedComp" "primeSieveSpeedComp. primesieveseq (10000000)" 10个循环,最好是3:530毫秒每循环

import primeSieveSpeedComp (primeSieveSpeedComp) “primeSieveSpeedComp.rwh_primes1(10000000)”10圈,3:494的最佳成绩 每回路Msec $ python -m timeit -s"import primeSieveSpeedComp" “primeSieveSpeedComp.rwh_primes2(10000000)”10圈,最好的3:375 每回路Msec

10^8

$ python -mtimeit -s"import primeSieveSpeedComp" "primeSieveSpeedComp. primesieveseq (100000000)" 10圈,最好的3:5.55秒每圈

import primeSieveSpeedComp (primeSieveSpeedComp) “primeSieveSpeedComp.rwh_primes1(100000000)”10圈,最好成绩是3:5.33 秒/循环 $ python -m timeit -s"import primeSieveSpeedComp" “primeSieveSpeedComp.rwh_primes2(100000000)”10圈,最好的3:3.95 秒/循环

10^9

$ python -mtimeit -s"import primeSieveSpeedComp" "primeSieveSpeedComp. primesieveseq (1000000000)" 10圈,最好的3圈:每圈61.2秒

$ python -mtimeit -n 3 -s"import primeSieveSpeedComp" “primeSieveSpeedComp.rwh_primes1(1000000000)”3圈,最好的3:97.8 秒/循环 $ python -m timeit -s"import primeSieveSpeedComp" “primeSieveSpeedComp.rwh_primes2(1000000000)”10个循环,3个最好: 每循环41.9秒

您可以将下面的代码复制到ubuntu primeSieveSpeedComp中以查看此测试。

def primeSieveSeq(MAX_Int):
    if MAX_Int > 5*10**8:
        import ctypes
        int16Array = ctypes.c_ushort * (MAX_Int >> 1)
        sieve = int16Array()
        #print 'uses ctypes "unsigned short int Array"'
    else:
        sieve = (MAX_Int >> 1) * [False]
        #print 'uses python list() of long long int'
    if MAX_Int < 10**8:
        sieve[4::3] = [True]*((MAX_Int - 8)/6+1)
        sieve[12::5] = [True]*((MAX_Int - 24)/10+1)
    r = [2, 3, 5]
    n = 0
    for i in xrange(int(MAX_Int**0.5)/30+1):
        n += 3
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 2
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 1
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 2
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 1
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 2
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 3
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
        n += 1
        if not sieve[n]:
            n2 = (n << 1) + 1
            r.append(n2)
            n2q = (n2**2) >> 1
            sieve[n2q::n2] = [True]*(((MAX_Int >> 1) - n2q - 1) / n2 + 1)
    if MAX_Int < 10**8:
        return [2, 3, 5]+[(p << 1) + 1 for p in [n for n in xrange(3, MAX_Int >> 1) if not sieve[n]]]
    n = n >> 1
    try:
        for i in xrange((MAX_Int-2*n)/30 + 1):
            n += 3
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 2
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 1
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 2
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 1
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 2
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 3
            if not sieve[n]:
                r.append((n << 1) + 1)
            n += 1
            if not sieve[n]:
                r.append((n << 1) + 1)
    except:
        pass
    return r

其他回答

下面是Eratosthenes的一个numpy版本,具有良好的复杂度(低于排序长度为n的数组)和向量化。与@unutbu相比,用46微秒就可以找到100万以下的所有质数。

import numpy as np 
def generate_primes(n):
    is_prime = np.ones(n+1,dtype=bool)
    is_prime[0:2] = False
    for i in range(int(n**0.5)+1):
        if is_prime[i]:
            is_prime[i**2::i]=False
    return np.where(is_prime)[0]

计时:

import time    
for i in range(2,10):
    timer =time.time()
    generate_primes(10**i)
    print('n = 10^',i,' time =', round(time.time()-timer,6))

>> n = 10^ 2  time = 5.6e-05
>> n = 10^ 3  time = 6.4e-05
>> n = 10^ 4  time = 0.000114
>> n = 10^ 5  time = 0.000593
>> n = 10^ 6  time = 0.00467
>> n = 10^ 7  time = 0.177758
>> n = 10^ 8  time = 1.701312
>> n = 10^ 9  time = 19.322478

如果你可以控制N,列出所有质数的最快方法就是预先计算它们。认真对待。预计算是一种被忽视的优化方法。

这里是最快的函数之一的两个更新版本(纯Python 3.6),

from itertools import compress

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

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

使用Sundaram的Sieve,我想我打破了pure-Python的记录:

def sundaram3(max_n):
    numbers = range(3, max_n+1, 2)
    half = (max_n)//2
    initial = 4

    for step in xrange(3, max_n+1, 2):
        for i in xrange(initial, half, step):
            numbers[i-1] = 0
        initial += 2*(step+1)

        if initial > half:
            return [2] + filter(None, numbers)

Comparasion:

C:\USERS>python -m timeit -n10 -s "import get_primes" "get_primes.get_primes_erat(1000000)"
10 loops, best of 3: 710 msec per loop

C:\USERS>python -m timeit -n10 -s "import get_primes" "get_primes.daniel_sieve_2(1000000)"
10 loops, best of 3: 435 msec per loop

C:\USERS>python -m timeit -n10 -s "import get_primes" "get_primes.sundaram3(1000000)"
10 loops, best of 3: 327 msec per loop

下面是一个使用python的列表推导式生成质数的有趣技术(但不是最有效的):

noprimes = [j for i in range(2, 8) for j in range(i*2, 50, i)]
primes = [x for x in range(2, 50) if x not in noprimes]