2

最近、素数ふるいを高速化するために、ここから bitarray モジュールをダウンロードしましたが、結果は悲惨です。

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

ふるいは、リストのbitarray2 倍以上遅いです。より良い配列の提案はありますか? 何でも python よりも高速である必要があるlistか、そう思いました。 最速ですが、時間がかかるのでnumpy.ones嫌いです。numpyimport

True私は基本的に、変更可能で、およびを保持できる高速なデータホルダーを探していFalseます。

4

1 に答える 1

2

実際のビットの設定とクリアは、bitarray では非常に高速です。作成しているのは、戻りリストの作成が遅いことです。範囲を反復してから各ビットをテストする代わりに、ビットを反復する bitarray のサポートを利用します。

これを試して：

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.47秒かかりますが、私のマシンでは約0.38秒で実行されます。

この質問を見たことがありますか？

私は 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. 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))

ビットの実際の設定/クリアについては、bitarray が gmpy2 よりも高速です。また、bitarray には、gmpy2 にはない多くの機能があります。ただし、bitarrary で、どのビットが設定またはクリアされているかのインデックスを取得するためのより高速な方法を見つけることができませんでした。

ところで、ベンチマーク関数 sieve2() および sieve3() は間違った結果を返します。あなたが行方不明if a[x]です。

于 2013-03-28T05:40:43.047 に答える