最快可变数据保持器

Question

我最近从这里下载了位数组模块，以获得更快的素筛，但是结果很糟糕。

from bitarray import bitarray
from numpy import ones
from timeit import timeit


def sieve1(n = 1000000):
    '''Sieve of Eratosthenes (bitarray)'''
    l = (n - 1)/2; a = bitarray(l); a.setall(True)
    for i in xrange(500):
        if a[i]:
            s = i+i+3; t = (s*s-3)/2; a[t:l:s] = False
    return [2] + [x+x+3 for x in xrange(l) if a[x]]

def sieve2(n = 1000000):
    '''Sieve of Eratosthenes (list)'''
    l = (n - 1)/2; a = [True] * l
    for i in xrange(500):
        if a[i]:
            s = i+i+3; t = (s*s-3)/2; u = l-t-1
            a[t:l:s] = [False] * (u/s + 1)
    return [2] + [x+x+3 for x in xrange(l)]

def sieve3(n = 1000000):
    '''Sieve of Eratosthenes (numpy.ones)'''    
    l = (n - 1)/2; a = ones(l, dtype=bool)
    for i in xrange(500):
        if a[i]:
            s = i+i+3; t = (s*s-3)/2; a[t:l:s] = False
    return [2] + [x+x+3 for x in xrange(l)]

print timeit(sieve1, number=10)
print timeit(sieve2, number=10)
print timeit(sieve3, number=10)

以下是结果-

1.59695601594
0.666230770593
0.523708537583

bitarray筛子的速度是列表的两倍多。有人对更好的数组有建议吗？任何事情都应该比python list更快，至少我是这么想的。numpy.ones是最快的，但我不喜欢numpy，因为import需要很长时间。

我基本上是在寻找一个快速的数据持有者，它是可变的，可以容纳True和False。

Answer 1

在位数组中，位的实际设置和清除非常快。正在使构建返回列表的速度变慢。与其遍历范围，然后测试每一位，不如利用位数组对遍历位的支持。

试试这个：

def bitarray_sieve(n = 1000000):
    '''Sieve of Eratosthenes (bitarray)'''
    l = (n - 1)//2; a = bitarray(l); a.setall(True)
    for i in range(500):
        if a[i]:
            s = i+i+3; t = (s*s-3)//2; a[t:l:s] = False
    return [2] + [x+x+3 for x,b in enumerate(a) if b]

它在我的机器上运行大约0.38秒，而列表版本大约需要0.47秒。

你看过这个问题吗？

我维护gmpy2，并增加了在整数中迭代位和设置/清除位的能力。

下面的示例大约需要0.16秒。

def gmpy2_sieve2(n=1000000):
    '''Sieve of Eratosthenes (gmpy2, version 2)'''
    l = (n - 1)//2; a = gmpy2.xbit_mask(l)
    for i in range(500):
        if a[i]:
            s = i+i+3; t = (s*s-3)//2; u = l-t-1
            a[t:l:s] = 0
    return [2] + [x+x+3 for x in a.iter_set()]

现在的瓶颈是计算x+x+3。下面的解决方案不会跳过筛分2，它需要两倍的内存，但它允许立即使用位位置。在我的机器上大约需要0.08秒：

def gmpy2_sieve(limit=1000000):
    '''Returns a generator that yields the prime numbers up to limit.

    Bits are set to 1 if their position is composite.'''

    sieve_limit = gmpy2.isqrt(limit) + 1
    limit += 1

    # Mark bit positions 0 and 1 as not prime.
    bitmap = gmpy2.xmpz(3)

    # Process 2 separately. This allows us to use p+p for the step size
    # when sieving the remaining primes.
    bitmap[4 : limit : 2] = -1

    # Sieve the remaining primes.
    for p in bitmap.iter_clear(3, sieve_limit):
        bitmap[p*p : limit : p+p] = -1

    return list(bitmap.iter_clear(2, limit))

对于位的实际设置/清除，位数组比gmpy2快。而且位数组具有gmpy2所缺乏的许多功能。但是，我无法在位元中找到一个更快的方法来获得设置或清除哪个位的索引。

顺便说一句，您的基准函数sieve2()和sieve3()返回不正确的结果；您缺少if a[x]。