Accélérer les opérations bitstring/bit en Python?
j'ai écrit un générateur de nombres premiers en utilisant tamis D'Eratosthène et Python 3.1. Le code fonctionne correctement et gracieusement à 0.32 secondes sur ideone.com pour générer des nombres premiers jusqu'à 1.000.000.
# from bitstring import BitString
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
flags = [False, False] + [True] * (limit - 2)
# flags = BitString(limit)
# Step through all the odd numbers
for i in range(3, limit, 2):
if flags[i] is False:
# if flags[i] is True:
continue
yield i
# Exclude further multiples of the current prime number
if i <= sub_limit:
for j in range(i*3, limit, i<<1):
flags[j] = False
# flags[j] = True
le problème est, je cours de la mémoire quand j'essaie de générer des nombres jusqu'à 1.000.000.000.
flags = [False, False] + [True] * (limit - 2)
MemoryError
comme vous pouvez l'imaginer, allouant 1 milliard de valeurs booléennes ( 1 byte 4 ou 8 bytes (voir commentaire) chacun en Python) est vraiment impossible, donc j'ai regardé dans bitstring . J'ai pensé qu'utiliser 1 bit pour chaque drapeau serait beaucoup plus efficace pour la mémoire. Cependant, la performance du programme a chuté drastiquement - 24 secondes d'exécution, pour le nombre premier jusqu'à 1.000.000. Ceci est probablement dû à l'implémentation interne de bitstring.
vous pouvez commenter / décommenter les trois lignes pour voir ce que j'ai changé pour utiliser BitString, comme l'extrait de code ci-dessus.
ma question Est, y a-t-il un moyen d'accélérer mon programme, avec ou sans bitstring?
Edit: s'il vous Plaît tester le code vous-même avant de poster. Je ne peux pas accepter des réponses qui fonctionnent plus lentement que mon code existant, naturellement.
Éditer à nouveau:
11 réponses
Il y a quelques petites optimisations pour votre version. En inversant les rôles de vrai et de faux, vous pouvez changer " if flags[i] is False: " en " if flags[i]: ". Et la valeur de départ pour la deuxième déclaration range peut être i*i au lieu de i*3 . Votre version originale prend 0.166 secondes sur mon système. Avec ces modifications, la version ci-dessous prend 0.156 secondes sur mon système.
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
flags = [True, True] + [False] * (limit - 2)
# Step through all the odd numbers
for i in range(3, limit, 2):
if flags[i]:
continue
yield i
# Exclude further multiples of the current prime number
if i <= sub_limit:
for j in range(i*i, limit, i<<1):
flags[j] = True
Ce n'aide pas votre problème de mémoire.
se déplaçant dans le monde des extensions C, j'ai utilisé la version de développement de gmpy . (Disclaimer: je suis l'un des responsables. La version de développement est appelée gmpy2 et supporte des entiers mutables appelés xmpz. En utilisant gmpy2 et le code suivant, j'ai un temps de fonctionnement de 0.140 secondes. Le temps de course pour une limite de 1.000.000.000 est de 158 secondes.
import gmpy2
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
# Actual number is 2*bit_position + 1.
oddnums = gmpy2.xmpz(1)
current = 0
while True:
current += 1
current = oddnums.bit_scan0(current)
prime = 2 * current + 1
if prime > limit:
break
yield prime
# Exclude further multiples of the current prime number
if prime <= sub_limit:
for j in range(2*current*(current+1), limit>>1, prime):
oddnums.bit_set(j)
poussant des optimisations, et sacrifiant la clarté, je reçois temps d'exécution de 0,107 et 123 secondes avec le code suivant:
import gmpy2
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
# Actual number is 2*bit_position + 1.
oddnums = gmpy2.xmpz(1)
f_set = oddnums.bit_set
f_scan0 = oddnums.bit_scan0
current = 0
while True:
current += 1
current = f_scan0(current)
prime = 2 * current + 1
if prime > limit:
break
yield prime
# Exclude further multiples of the current prime number
if prime <= sub_limit:
list(map(f_set,range(2*current*(current+1), limit>>1, prime)))
modifier: basé sur cet exercice, j'ai modifié gmpy2 pour accepter xmpz.bit_set(iterator) . En utilisant le code suivant, le temps d'exécution pour tous les nombres premiers moins 1.000.000.000 est de 56 secondes pour Python 2.7 et 74 secondes pour Python 3.2. (Tel que noté dans les commentaires, xrange est plus rapide que range .)
import gmpy2
try:
range = xrange
except NameError:
pass
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
oddnums = gmpy2.xmpz(1)
f_scan0 = oddnums.bit_scan0
current = 0
while True:
current += 1
current = f_scan0(current)
prime = 2 * current + 1
if prime > limit:
break
yield prime
if prime <= sub_limit:
oddnums.bit_set(iter(range(2*current*(current+1), limit>>1, prime)))
Edit #2: encore Un essai! J'ai modifié gmpy2 pour accepter xmpz.bit_set(slice) . En utilisant le code suivant, le temps d'exécution pour tous les nombres premiers moins 1.000.000.000 est d'environ 40 secondes pour Python 2.7 et Python 3.2.
from __future__ import print_function
import time
import gmpy2
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
flags = gmpy2.xmpz(1)
# pre-allocate the total length
flags.bit_set((limit>>1)+1)
f_scan0 = flags.bit_scan0
current = 0
while True:
current += 1
current = f_scan0(current)
prime = 2 * current + 1
if prime > limit:
break
yield prime
if prime <= sub_limit:
flags.bit_set(slice(2*current*(current+1), limit>>1, prime))
start = time.time()
result = list(prime_numbers(1000000000))
print(time.time() - start)
Edit #3: j'ai mis à jour gmpy2 pour supporter correctement le tranchage au niveau des bits d'une xmpz. Aucun changement dans la performance, mais une API très agréable. J'ai fait quelques ajustements et j'ai eu le temps d'environ 37 secondes. (Voir Modifier # 4 pour les changements dans gmpy2 2.0.0b1.)
from __future__ import print_function
import time
import gmpy2
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
sub_limit = int(limit**0.5)
flags = gmpy2.xmpz(1)
flags[(limit>>1)+1] = True
f_scan0 = flags.bit_scan0
current = 0
prime = 2
while prime <= sub_limit:
yield prime
current += 1
current = f_scan0(current)
prime = 2 * current + 1
flags[2*current*(current+1):limit>>1:prime] = True
while prime <= limit:
yield prime
current += 1
current = f_scan0(current)
prime = 2 * current + 1
start = time.time()
result = list(prime_numbers(1000000000))
print(time.time() - start)
Edit #4: j'ai apporté quelques modifications en gmpy2 2.0.0b1 qui casse l'exemple précédent. gmpy2 ne traite plus Vrai qu'une valeur spéciale qui fournit une source infinie de 1 bits. -1 doit être utilisé à la place.
from __future__ import print_function
import time
import gmpy2
def prime_numbers(limit=1000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
sub_limit = int(limit**0.5)
flags = gmpy2.xmpz(1)
flags[(limit>>1)+1] = 1
f_scan0 = flags.bit_scan0
current = 0
prime = 2
while prime <= sub_limit:
yield prime
current += 1
current = f_scan0(current)
prime = 2 * current + 1
flags[2*current*(current+1):limit>>1:prime] = -1
while prime <= limit:
yield prime
current += 1
current = f_scan0(current)
prime = 2 * current + 1
start = time.time()
result = list(prime_numbers(1000000000))
print(time.time() - start)
Edit #5: j'ai apporté quelques améliorations à gmpy2 2.0.0b2. Vous pouvez maintenant itérer tous les bits qui sont définis ou clairs. Le temps de fonctionnement s'est amélioré d'environ 30%.
from __future__ import print_function
import time
import gmpy2
def sieve(limit=1000000):
'''Returns a generator that yields the prime numbers up to limit.'''
# Increment by 1 to account for the fact that slices do not include
# the last index value but we do want to include the last value for
# calculating a list of primes.
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 bitmap.iter_clear(2, limit)
if __name__ == "__main__":
start = time.time()
result = list(sieve(1000000000))
print(time.time() - start)
print(len(result))
ok, voici un benchmarking complet que j'ai fait ce soir pour voir quel code fonctionne le plus vite. Espérons que quelqu'un trouvera cette liste utile. J'ai omis tout ce qui prend plus de 30 secondes à remplir sur ma machine.
je tiens à remercier tous ceux qui ont contribué. J'ai appris beaucoup de choses grâce à vos efforts, et j'espère que vous aussi.
ma machine: AMD ZM-86, 2.40 Ghz Dual-Core, avec 4 Go de RAM. C'est un portable HP Touchsmart Tx2. Notez que même si j'ai peut-être fait un lien avec un pastebin, , j'ai comparé tout ce qui suit sur ma propre machine.
j'ajouterai le benchmark gmpy2 une fois que je serai capable de le construire.
tous les benchmarks sont testés en python 2.6 x86
retour d'une liste de nombres premiers n jusqu'à 1.000.000: ( utilisant Python générateurs)
de Sebastian numpy générateur de version (mise à jour) - 121 ms @
tamis de Mark + roue - 154 ms
de Robert version avec découpage 1519210920" - 159 ms
Ma version améliorée avec de découpage - 205 ms
Numpy générateur avec énumérer - 249 ms @
tamis de base de Mark - 317 ms
amélioration de casevh sur mon original solution - 343 ms
Mon modifiée numpy générateur de solution - 407 ms
Ma méthode originale dans le question - 409 ms
Bitarray Solution - 414 ms @
Python pur avec bytearray - 1394 ms @
Scott's BitString solution - 6659 ms @
'@' signifie que cette méthode est capable de générer jusqu'à n < 1 000 000 000 d', sur mon installation de machine.
en outre, si vous n'avez pas besoin de génératrice et je veux juste toute la liste tout de suite:
numpy solution from RosettaCode - 32 ms @
(la solution de numpy est également capable de générer jusqu'à 1 milliard, ce qui a pris 61.6259 secondes. Je pense que la mémoire a été échangée une fois, d'où le double temps.)
OK, donc c'est ma deuxième réponse, mais comme la vitesse est essentielle, j'ai pensé que je devais mentionner le module bitarray - même si c'est bitstring nemesis :). Il est idéalement adapté à ce cas car non seulement c'est une extension C (et donc plus rapide que Python pur a un espoir d'être), mais il supporte aussi les assignations de tranches. Il n'est pas encore disponible pour Python 3.
Je n'ai même pas essayé d'optimiser cela, je j'ai réécrit la version bitstring. Sur ma machine, j'ai 0,16 seconde pour des nombres premiers inférieurs à un million.
pour un milliard, il fonctionne parfaitement et se termine en 2 minutes et 31 secondes.
import bitarray
def prime_bitarray(limit=1000000):
yield 2
flags = bitarray.bitarray(limit)
flags.setall(False)
sub_limit = int(limit**0.5)
for i in xrange(3, limit, 2):
if not flags[i]:
yield i
if i <= sub_limit:
flags[3*i:limit:i*2] = True
question connexe.
la manière la plus rapide de lister tous les nombres premiers ci-dessous N en python
Salut, je suis trop à la recherche d'un code en python pour générer des nombres premiers jusqu'à 10 9 aussi vite que je peux, ce qui est difficile à cause du problème de mémoire. jusqu'à présent, j'ai créé ceci pour générer des nombres premiers jusqu'à 10 6 & 10*7 (0,171 s et 1,764 s respectivement dans mon ancienne machine), qui semble être un peu plus rapide(au moins dans ma ordinateur) que "ma version améliorée avec tranchage" du post précédent, probablement parce que j'utilise n//i-i +1 (Voir ci-dessous)à la place de Len analogue(drapeaux[i2::i<<1]) dans l'autre code. s'il vous plaît corrigez-moi si je me trompe. Toute suggestion d'amélioration sont les bienvenus.
def primes(n):
""" Returns a list of primes < n """
sieve = [True] * n
for i in xrange(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 xrange(3,n,2) if sieve[i]]
dans L'un de ses codes Xavier utilise les drapeaux[i2::i<<1] et len(drapeaux[i2::i<<1]), j'ai calculé la taille explicite en utilisant le fait qu'entre i I..n Nous avons (N-I I)//2*I multiples de 2*I parce que nous voulons Compter i i aussi nous sommons 1 ainsi Len(tamis[i i::2*i]) égal (n-i*i)//(2*i) +1 cela rend le code plus rapide. Un générateur de base basé sur le code ci-dessus serait:
def primesgen(n):
""" Generates all primes <= n """
sieve = [True] * n
yield 2
for i in xrange(3,int(n**0.5)+1,2):
if sieve[i]:
yield i
sieve[i*i::2*i] = [False]*((n-i*i-1)/(2*i)+1)
for i in xrange(i+2,n,2):
if sieve[i]: yield i
avec un peu de modification on peut écrire une version légèrement plus lente du code ci-dessus qui commence avec un tamis moitié de la taille tamis = [True] * (n//2) et fonctionne pour le même N. Je ne sais pas comment vous réfléchissez dans la question de la mémoire. Comme un exemple de mise en œuvre ici est la version modifiée du code num PY rosetta (peut-être plus rapide) à partir d'un tamis de la moitié de la taille.
import numpy
def primesfrom3to(n):
""" Returns a array of primes, 3 <= p < n """
sieve = numpy.ones(n/2, dtype=numpy.bool)
for i in xrange(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
un générateur plus rapide et plus économe en mémoire serait:
import numpy
def primesgen1(n):
""" Input n>=6, Generates all primes < n """
sieve1 = numpy.ones(n/6+1, dtype=numpy.bool)
sieve5 = numpy.ones(n/6 , dtype=numpy.bool)
sieve1[0] = False
yield 2; yield 3;
for i in xrange(int(n**0.5)/6+1):
if sieve1[i]:
k=6*i+1; yield k;
sieve1[ ((k*k)/6) ::k] = False
sieve5[(k*k+4*k)/6::k] = False
if sieve5[i]:
k=6*i+5; yield k;
sieve1[ ((k*k)/6) ::k] = False
sieve5[(k*k+2*k)/6::k] = False
for i in xrange(i+1,n/6):
if sieve1[i]: yield 6*i+1
if sieve5[i]: yield 6*i+5
if n%6 > 1:
if sieve1[i+1]: yield 6*(i+1)+1
ou avec un peu plus de code:
import numpy
def primesgen(n):
""" Input n>=30, Generates all primes < n """
size = n/30 + 1
sieve01 = numpy.ones(size, dtype=numpy.bool)
sieve07 = numpy.ones(size, dtype=numpy.bool)
sieve11 = numpy.ones(size, dtype=numpy.bool)
sieve13 = numpy.ones(size, dtype=numpy.bool)
sieve17 = numpy.ones(size, dtype=numpy.bool)
sieve19 = numpy.ones(size, dtype=numpy.bool)
sieve23 = numpy.ones(size, dtype=numpy.bool)
sieve29 = numpy.ones(size, dtype=numpy.bool)
sieve01[0] = False
yield 2; yield 3; yield 5;
for i in xrange(int(n**0.5)/30+1):
if sieve01[i]:
k=30*i+1; yield k;
sieve01[ (k*k)/30::k] = False
sieve07[(k*k+ 6*k)/30::k] = False
sieve11[(k*k+10*k)/30::k] = False
sieve13[(k*k+12*k)/30::k] = False
sieve17[(k*k+16*k)/30::k] = False
sieve19[(k*k+18*k)/30::k] = False
sieve23[(k*k+22*k)/30::k] = False
sieve29[(k*k+28*k)/30::k] = False
if sieve07[i]:
k=30*i+7; yield k;
sieve01[(k*k+ 6*k)/30::k] = False
sieve07[(k*k+24*k)/30::k] = False
sieve11[(k*k+16*k)/30::k] = False
sieve13[(k*k+12*k)/30::k] = False
sieve17[(k*k+ 4*k)/30::k] = False
sieve19[ (k*k)/30::k] = False
sieve23[(k*k+22*k)/30::k] = False
sieve29[(k*k+10*k)/30::k] = False
if sieve11[i]:
k=30*i+11; yield k;
sieve01[ (k*k)/30::k] = False
sieve07[(k*k+ 6*k)/30::k] = False
sieve11[(k*k+20*k)/30::k] = False
sieve13[(k*k+12*k)/30::k] = False
sieve17[(k*k+26*k)/30::k] = False
sieve19[(k*k+18*k)/30::k] = False
sieve23[(k*k+ 2*k)/30::k] = False
sieve29[(k*k+ 8*k)/30::k] = False
if sieve13[i]:
k=30*i+13; yield k;
sieve01[(k*k+24*k)/30::k] = False
sieve07[(k*k+ 6*k)/30::k] = False
sieve11[(k*k+ 4*k)/30::k] = False
sieve13[(k*k+18*k)/30::k] = False
sieve17[(k*k+16*k)/30::k] = False
sieve19[ (k*k)/30::k] = False
sieve23[(k*k+28*k)/30::k] = False
sieve29[(k*k+10*k)/30::k] = False
if sieve17[i]:
k=30*i+17; yield k;
sieve01[(k*k+ 6*k)/30::k] = False
sieve07[(k*k+24*k)/30::k] = False
sieve11[(k*k+26*k)/30::k] = False
sieve13[(k*k+12*k)/30::k] = False
sieve17[(k*k+14*k)/30::k] = False
sieve19[ (k*k)/30::k] = False
sieve23[(k*k+ 2*k)/30::k] = False
sieve29[(k*k+20*k)/30::k] = False
if sieve19[i]:
k=30*i+19; yield k;
sieve01[ (k*k)/30::k] = False
sieve07[(k*k+24*k)/30::k] = False
sieve11[(k*k+10*k)/30::k] = False
sieve13[(k*k+18*k)/30::k] = False
sieve17[(k*k+ 4*k)/30::k] = False
sieve19[(k*k+12*k)/30::k] = False
sieve23[(k*k+28*k)/30::k] = False
sieve29[(k*k+22*k)/30::k] = False
if sieve23[i]:
k=30*i+23; yield k;
sieve01[(k*k+24*k)/30::k] = False
sieve07[(k*k+ 6*k)/30::k] = False
sieve11[(k*k+14*k)/30::k] = False
sieve13[(k*k+18*k)/30::k] = False
sieve17[(k*k+26*k)/30::k] = False
sieve19[ (k*k)/30::k] = False
sieve23[(k*k+ 8*k)/30::k] = False
sieve29[(k*k+20*k)/30::k] = False
if sieve29[i]:
k=30*i+29; yield k;
sieve01[ (k*k)/30::k] = False
sieve07[(k*k+24*k)/30::k] = False
sieve11[(k*k+20*k)/30::k] = False
sieve13[(k*k+18*k)/30::k] = False
sieve17[(k*k+14*k)/30::k] = False
sieve19[(k*k+12*k)/30::k] = False
sieve23[(k*k+ 8*k)/30::k] = False
sieve29[(k*k+ 2*k)/30::k] = False
for i in xrange(i+1,n/30):
if sieve01[i]: yield 30*i+1
if sieve07[i]: yield 30*i+7
if sieve11[i]: yield 30*i+11
if sieve13[i]: yield 30*i+13
if sieve17[i]: yield 30*i+17
if sieve19[i]: yield 30*i+19
if sieve23[i]: yield 30*i+23
if sieve29[i]: yield 30*i+29
i += 1
mod30 = n%30
if mod30 > 1:
if sieve01[i]: yield 30*i+1
if mod30 > 7:
if sieve07[i]: yield 30*i+7
if mod30 > 11:
if sieve11[i]: yield 30*i+11
if mod30 > 13:
if sieve13[i]: yield 30*i+13
if mod30 > 17:
if sieve17[i]: yield 30*i+17
if mod30 > 19:
if sieve19[i]: yield 30*i+19
if mod30 > 23:
if sieve23[i]: yield 30*i+23
Ps: Si vous avez des suggestions sur la façon d'accélérer ce code, n'importe quoi de l'ordre de changement des opérations au pré-calcul n'importe quoi, s'il vous plaît commenter.
Voici une version que j'ai écrite il y a un certain temps; il pourrait être intéressant de comparer avec la vôtre pour la vitesse. Il ne fait rien au sujet des problèmes d'espace, cependant (en fait, ils sont probablement pires qu'avec votre version).
from math import sqrt
def basicSieve(n):
"""Given a positive integer n, generate the primes < n."""
s = [1]*n
for p in xrange(2, 1+int(sqrt(n-1))):
if s[p]:
a = p*p
s[a::p] = [0] * -((a-n)//p)
for p in xrange(2, n):
if s[p]:
yield p
j'ai des versions plus rapides, utilisant une roue, mais elles sont beaucoup plus compliquées. La source originale est ici .
ok, voici la version avec une roue. wheelSieve est la principale entrée point.
from math import sqrt
from bisect import bisect_left
def basicSieve(n):
"""Given a positive integer n, generate the primes < n."""
s = [1]*n
for p in xrange(2, 1+int(sqrt(n-1))):
if s[p]:
a = p*p
s[a::p] = [0] * -((a-n)//p)
for p in xrange(2, n):
if s[p]:
yield p
class Wheel(object):
"""Class representing a wheel.
Attributes:
primelimit -> wheel covers primes < primelimit.
For example, given a primelimit of 6
the wheel primes are 2, 3, and 5.
primes -> list of primes less than primelimit
size -> product of the primes in primes; the modulus of the wheel
units -> list of units modulo size
phi -> number of units
"""
def __init__(self, primelimit):
self.primelimit = primelimit
self.primes = list(basicSieve(primelimit))
# compute the size of the wheel
size = 1
for p in self.primes:
size *= p
self.size = size
# find the units by sieving
units = [1] * self.size
for p in self.primes:
units[::p] = [0]*(self.size//p)
self.units = [i for i in xrange(self.size) if units[i]]
# number of units
self.phi = len(self.units)
def to_index(self, n):
"""Compute alpha(n), where alpha is an order preserving map
from the set of units modulo self.size to the nonnegative integers.
If n is not a unit, the index of the first unit greater than n
is given."""
return bisect_left(self.units, n % self.size) + n // self.size * self.phi
def from_index(self, i):
"""Inverse of to_index."""
return self.units[i % self.phi] + i // self.phi * self.size
def wheelSieveInner(n, wheel):
"""Given a positive integer n and a wheel, find the wheel indices of
all primes that are less than n, and that are units modulo the
wheel modulus.
"""
# renaming to avoid the overhead of attribute lookups
U = wheel.units
wS = wheel.size
# inverse of unit map
UI = dict((u, i) for i, u in enumerate(U))
nU = len(wheel.units)
sqroot = 1+int(sqrt(n-1)) # ceiling of square root of n
# corresponding index (index of next unit up)
sqrti = bisect_left(U, sqroot % wS) + sqroot//wS*nU
upper = bisect_left(U, n % wS) + n//wS*nU
ind2 = bisect_left(U, 2 % wS) + 2//wS*nU
s = [1]*upper
for i in xrange(ind2, sqrti):
if s[i]:
q = i//nU
u = U[i%nU]
p = q*wS+u
u2 = u*u
aq, au = (p+u)*q+u2//wS, u2%wS
wp = p * nU
for v in U:
# eliminate entries corresponding to integers
# congruent to p*v modulo p*wS
uvr = u*v%wS
m = aq + (au > uvr)
bot = (m + (q*v + u*v//wS - m) % p) * nU + UI[uvr]
s[bot::wp] = [0]*-((bot-upper)//wp)
return s
def wheelSieve(n, wheel=Wheel(10)):
"""Given a positive integer n, generate the list of primes <= n."""
n += 1
wS = wheel.size
U = wheel.units
nU = len(wheel.units)
s = wheelSieveInner(n, wheel)
return (wheel.primes[:bisect_left(wheel.primes, n)] +
[p//nU*wS + U[p%nU] for p in xrange(bisect_left(U, 2 % wS)
+ 2//wS*nU, len(s)) if s[p]])
ce Que la roue est: bien, vous savez que (hormis 2), tous les nombres premiers sont impairs, donc la plupart des tamis manquer tous les nombres pairs. De même, vous pouvez aller un peu plus loin et noter que tous les nombres premiers (sauf 2 et 3) sont congruents à 1 ou 5 modulo 6 (== 2 * 3), de sorte que vous pouvez vous en tirer avec seulement le stockage des entrées pour ces nombres dans votre tamis. L'étape suivante consiste à noter que tous les nombres premiers (sauf 2, 3 et 5) sont congruents à l'un des nombres suivants: 1, 7, 11, 13, 17, 19, 23, 29 (modulo 30) (ici 30 == 2*3*5), et ainsi de suite.
une amélioration de la vitesse que vous pouvez faire en utilisant bitstring est d'utiliser la méthode 'set' lorsque vous excluez des multiples du nombre courant.
ainsi la section vitale devient
for i in range(3, limit, 2):
if flags[i]:
yield i
if i <= sub_limit:
flags.set(1, range(i*3, limit, i*2))
sur ma machine cela fonctionne environ 3 fois plus vite que l'original.
mon autre pensée était d'utiliser la chaîne de bits pour représenter seulement les nombres impairs. Vous pouvez alors trouver les bits non connectés en utilisant flags.findall([0]) qui renvoie un générateur. Pas sûr si ce serait beaucoup plus rapide (Trouver des modèles de bits n'est pas terriblement facile à faire efficacement).
[divulgation complète: j'ai écrit le module bitstring, donc j'ai une certaine fierté en jeu ici!]
en comparaison, j'ai aussi retiré les intestins de la méthode bitstring de sorte qu'il le fait de la même manière, mais sans vérification de la portée, etc. Je pense que cela donne une limite inférieure raisonnable pour Python pur qui fonctionne pour un milliard d'éléments (sans changer l'algorithme, ce qui, je pense, c'est de la triche!)
def prime_pure(limit=1000000):
yield 2
flags = bytearray((limit + 7) // 8)
sub_limit = int(limit**0.5)
for i in xrange(3, limit, 2):
byte, bit = divmod(i, 8)
if not flags[byte] & (128 >> bit):
yield i
if i <= sub_limit:
for j in xrange(i*3, limit, i*2):
byte, bit = divmod(j, 8)
flags[byte] |= (128 >> bit)
sur ma machine, cela prend environ 0,62 seconde pour un million d'éléments, ce qui signifie que c'est environ un quart de la vitesse de la réponse bitarray. Je pense que c'est tout à fait raisonnable pour Python pur.
les types entiers de Python peuvent être de taille arbitraire, donc vous ne devriez pas avoir besoin d'une bibliothèque de bitstring intelligente pour représenter un ensemble de bits, juste un seul nombre.
voici le code. Il gère une limite de 1.000.000 avec facilité, et peut même gérer 10.000.000 sans se plaindre trop:
def multiples_of(n, step, limit):
bits = 1 << n
old_bits = bits
max = 1 << limit
while old_bits < max:
old_bits = bits
bits += bits << step
step *= 2
return old_bits
def prime_numbers(limit=10000000):
'''Prime number generator. Yields the series
2, 3, 5, 7, 11, 13, 17, 19, 23, 29 ...
using Sieve of Eratosthenes.
'''
yield 2
sub_limit = int(limit**0.5)
flags = ((1 << (limit - 2)) - 1) << 2
# Step through all the odd numbers
for i in xrange(3, limit, 2):
if not (flags & (1 << i)):
continue
yield i
# Exclude further multiples of the current prime number
if i <= sub_limit:
flags &= ~multiples_of(i * 3, i << 1, limit)
la fonction multiples_of évite le coût de réglage de chaque multiple individuellement. Il règle le bit initial, le déplace par l'étape initiale ( i << 1 ) et ajoute le résultat à lui-même. Il Double alors le pas, de sorte que le prochain décalage copie les deux bits, et ainsi de suite jusqu'à ce qu'il atteigne la limite. Cela règle Tous les multiples d'un nombre jusqu'à la limite dans le temps o(log(limite)).
l'une des options que vous voudrez peut-être Regarder est de compiler un module C/C++ pour avoir un accès direct aux fonctionnalités de "bit-twiddling". AFAIK vous pouvez écrire quelque chose de cette nature et ne retourner les résultats qu'à la fin des calculs effectués en C/C++. Maintenant que je tape ceci, vous pouvez aussi regarder numpy qui stocke des tableaux en tant que c natif pour la vitesse. Je ne sais pas si ce sera plus rapide que le module bitstring!
une autre option vraiment stupide, mais cela peut être d'une aide si vous avez vraiment besoin d'une grande liste de nombres premiers disponibles très rapidement. Par exemple, si vous avez besoin d'eux comme un outil pour répondre aux problèmes du projet Euler (si le problème est juste l'optimisation de votre code comme un jeu d'esprit, il est hors de propos).
utilisez n'importe quelle solution lente pour générer la liste et la sauvegarder dans un fichier source python, dit primes.py qui contiendrait juste:
primes = [ a list of a million primes numbers here ]
maintenant pour les utiliser vous venez vous devez faire import primes et vous les obtenez avec une empreinte mémoire minimale à la vitesse du disque IO. La complexité est le nombre de nombres premiers: -)
même si vous avez utilisé une solution mal optimisée pour générer cette liste, cela ne sera fait qu'une seule fois et cela n'a pas beaucoup d'importance.
vous pourriez probablement le rendre encore plus rapide en utilisant pickle / unpickle, mais vous avez l'idée...
vous pourriez utiliser un tamis segmenté de Eratosthènes. La mémoire utilisée pour chaque segment est réutilisée pour le segment suivant.
voici un code Python3 qui utilise moins de mémoire que les solutions bitarray/bitstring postées à ce jour et environ 1/8 de la mémoire de primesgen () de Robert William Hanks, tout en courant plus vite que primesgen() (légèrement plus rapide à 1.000.000, en utilisant 37KB de mémoire, 3x plus rapide que primesgen() à 1.000.000.000 en utilisant moins de 34MB). Augmenter la taille du morceau (variable chunk dans le code) utilise plus de mémoire mais accélère le programme, jusqu'à une limite -- j'ai choisi une valeur de sorte que sa contribution à la mémoire est sous environ 10% de n//30 octets de tamis. Il n'est pas aussi efficace en mémoire que le générateur infini de Will Ness ( https://stackoverflow.com/a/19391111/5439078 ) parce qu'il enregistre (et retourne à la fin, sous forme comprimée) tous les nombres premiers calculés.
cela devrait fonctionner correctement tant que le calcul à la racine carrée est exact (environ 2**51 si Python utilise des doubles 64 bits). Cependant, vous ne devriez pas utiliser ce programme pour trouver des nombres premiers aussi grands!
(Je n'ai pas recalculé les offsets, je les ai juste pris du code de Robert William Hanks. Merci Robert!)
import numpy as np
def prime_30_gen(n):
""" Input n, int-like. Yields all primes < n as Python ints"""
wheel = np.array([2,3,5])
yield from wheel[wheel < n].tolist()
powers = 1 << np.arange(8, dtype='u1')
odds = np.array([1, 7, 11, 13, 17, 19, 23, 29], dtype='i8')
offsets=np.array([[0,6,10,12,16,18,22,28],[6,24,16,12,4,0,22,10],
[0,6,20,12,26,18,2,8], [24,6,4,18,16,0,28,10],
[6,24,26,12,14,0,2,20], [0,24,10,18,4,12,28,22],
[24,6,14,18,26,0,8,20], [0,24,20,18,14,12,8,2]],
dtype='i8')
offsets = {d:f for d,f in zip(odds, offsets)}
sieve = 255 * np.ones((n + 28) // 30, dtype='u1')
if n % 30 > 1: sieve[-1] >>= 8 - odds.searchsorted(n % 30)
sieve[0] &= ~1
for i in range((int((n - 1)**0.5) + 29) // 30):
for odd in odds[(sieve[i] & powers).nonzero()[0]]:
k = i * 30 + odd
yield int(k)
for clear_bit, off in zip(~powers, offsets[odd]):
sieve[(k * (k + off)) // 30 :: k] &= clear_bit
chunk = min(1 + (n >> 13), 1<<15)
for i in range(i + 1, len(sieve), chunk):
a = (sieve[i : i + chunk, np.newaxis] & powers)
a = np.flatnonzero(a.astype('bool')) + i * 8
a = (a // 8 * 30 + odds[a & 7]).tolist()
yield from a
return sieve
Side note: si vous voulez une vitesse réelle et que vous avez les 2 Go de RAM nécessaires pour la liste des nombres premiers à 10* * 9, alors vous devez utiliser typrimesieve (sur https://pypi.python.org / , en utilisant primesieve http://primesieve.org / ).