while - numeros primos recursividad python

¿Cómo implementar un generador infinito eficiente de números primos en Python?

Question

Esto no es una tarea, solo tengo curiosidad.

INFINITE es la palabra clave aquí.

Deseo usarlo como para p en primos (). Creo que esta es una función incorporada en Haskell.

Entonces, la respuesta no puede ser tan ingenua como "Just do a Sieve".

En primer lugar, no sabe cuántos primos consecutivos se consumirán. Bueno, supongamos que puedes inventar 100 de ellos a la vez. ¿Usarías el mismo enfoque Sieve así como la frecuencia de la fórmula de los números primos?

Prefiero el enfoque no concurrente.

Gracias por leer (y escribir;))!

Answer 1

"Si he visto más ..."

La función erat2 del libro de cocina se puede acelerar aún más (en aproximadamente un 20-25%):

erat2a

import itertools as it
def erat2a( ):
    D = {  }
    yield 2
    for q in it.islice(it.count(3), 0, None, 2):
        p = D.pop(q, None)
        if p is None:
            D[q*q] = q
            yield q
        else:
            # old code here:
            # x = p + q
            # while x in D or not (x&1):
            #     x += p
            # changed into:
            x = q + 2*p
            while x in D:
                x += 2*p
            D[x] = p

La comprobación not (x&1) verifica que x es impar. Sin embargo, dado que tanto q como p son impares, al sumar 2*p se evitan la mitad de los pasos junto con la prueba de imparidad.

erat3

Si a uno no le importa un poco de extravagancia extra, erat2 se puede acelerar un 35-40% con los siguientes cambios (NB: necesita Python 2.7+ o Python 3+ debido a la función de itertools.compress ):

import itertools as it
def erat3( ):
    D = { 9: 3, 25: 5 }
    yield 2
    yield 3
    yield 5
    MASK= 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0,
    MODULOS= frozenset( (1, 7, 11, 13, 17, 19, 23, 29) )

    for q in it.compress(
            it.islice(it.count(7), 0, None, 2),
            it.cycle(MASK)):
        p = D.pop(q, None)
        if p is None:
            D[q*q] = q
            yield q
        else:
            x = q + 2*p
            while x in D or (x%30) not in MODULOS:
                x += 2*p
            D[x] = p

La función erat3 aprovecha el hecho de que todas las primos (a excepción de 2, 3, 5) módulo 30 resultan en solo ocho números: los incluidos en los MODULOS frozenset. Por lo tanto, después de ceder los tres primos iniciales, comenzamos desde 7 y trabajamos solo con los candidatos.
El filtrado candidato usa la función itertools.compress ; la "magia" está en la secuencia MASK ; MASK tiene 15 elementos (hay 15 números impares en cada 30 números, elegidos por la función itertools.islice ) con un 1 para cada posible candidato, comenzando desde 7. El ciclo se repite tal como lo especifica la función itertools.cycle .
La introducción del filtrado candidato necesita otra modificación: el or (x%30) not in MODULOS . El algoritmo erat2 procesó todos los números impares; ahora que el algoritmo erat3 procesa solo candidatos r30, debemos asegurarnos de que todas las D.keys() solo puedan ser candidatos falsos.

Puntos de referencia

Resultados

En un servidor Atom 330 Ubuntu 9.10, versiones 2.6.4 y 3.1.1+:

$ testit
up to 8192
==== python2 erat2 ====
100 loops, best of 3: 18.6 msec per loop
==== python2 erat2a ====
100 loops, best of 3: 14.5 msec per loop
==== python2 erat3 ====
Traceback (most recent call last):
…
AttributeError: ''module'' object has no attribute ''compress''
==== python3 erat2 ====
100 loops, best of 3: 19.2 msec per loop
==== python3 erat2a ====
100 loops, best of 3: 14.1 msec per loop
==== python3 erat3 ====
100 loops, best of 3: 11.7 msec per loop

En un servidor doméstico AMD Geode LX Gentoo, Python 2.6.5 y 3.1.2:

$ testit
up to 8192
==== python2 erat2 ====
10 loops, best of 3: 104 msec per loop
==== python2 erat2a ====
10 loops, best of 3: 81 msec per loop
==== python2 erat3 ====
Traceback (most recent call last):
…
AttributeError: ''module'' object has no attribute ''compress''
==== python3 erat2 ====
10 loops, best of 3: 116 msec per loop
==== python3 erat2a ====
10 loops, best of 3: 82 msec per loop
==== python3 erat3 ====
10 loops, best of 3: 66 msec per loop

Código de referencia

Un módulo primegen.py contiene las erat2 , erat2a y erat3 . Aquí sigue el guión de prueba:

#!/bin/sh
max_num=${1:-8192}
echo up to $max_num
for python_version in python2 python3
do
    for function in erat2 erat2a erat3
    do
        echo "==== $python_version $function ===="
        $python_version -O -m timeit -c /
        -s  "import itertools as it, functools as ft, operator as op, primegen; cmp= ft.partial(op.ge, $max_num)" /
            "next(it.dropwhile(cmp, primegen.$function()))"
    done
done

Answer 2

Aquí hay un generador infinito bastante rápido, escrito en Python2 pero que se ajusta fácilmente a Python3. Para usarlo para agregar los números primos hasta 10 ** 9, use lo siguiente:

from itertools import takewhile
from functools import partial
from operator import gt
print (sum(takewhile(partial(gt, 10**9), prime_gen_inf())))

Es un tamiz segmentado, más rápido pero obviamente menos elegante que el algoritmo de Will Ness.

from operator import mul
from functools import reduce
def prod(x): return reduce(mul, x, 1)


def build_sieve(wheel):
    w = prod(wheel)
    w_phi = prod([p-1 for p in wheel])
    rems = [a for a in range(w) if all(a % p for p in wheel)]
    assert len(rems) == w_phi
    inv = {a:pow(a, w_phi - 1, w) for a in rems}
    try:
        known_p = wheel + rems[1 : rems.index(rems[1]*rems[1])]
    except ValueError:
        known_p = wheel + rems[1:]
    return wheel, w, w_phi, rems, inv, known_p

#Adjust the chunk variable based on your computer''s architecture.
#
#Adjust the line with #! if you don''t need "true" infinite.  If you don''t need
#primes larger than 1<<32, use array(''H'', []), if 1<<64 use ''L'', if 1<<128 (in
#Python3) use ''Q'', otherwise use empty list [].
#To save memory, comment out the lines with #*, and uncomment the commented-out
#lines 
import itertools
from itertools import islice, count, compress, izip
chain_f = itertools.chain.from_iterable
from array import array
def prime_gen_inf(chunk=250000, sieve_info = build_sieve([2,3,5,7])):
    """    Indefinitely yields primes    """
    wheel, w, w_phi, rems, inv, known_p = sieve_info
    for p in known_p: yield p
    new_n = 0;
    while True:
        size = min(chunk, (p * p - new_n) / w)
        sieve = bytearray([1]) * size * w_phi
        n, new_n = new_n, new_n + size * w
        if not n:
            zero = bytearray([0])
            seen = len(known_p) - len(wheel) + 1
            sieve[:seen:1] = zero * seen
            p_gen = islice(prime_gen_inf(), len(wheel), None)
            new_p = next(p_gen)
            ps = []                                         #! array(''H'', [])
            p_invs = bytearray([])                                         #*
        while new_p * new_p < new_n:
            ps.append(new_p)
            p_invs.append(inv[new_p % w])                                  #*
            new_p = next(p_gen)
        for p, p_inv, modp in izip(ps, p_invs, [-n % p for p in ps]):      #*
            s = [(modp + p * (p_inv * (r - modp) % w)) / w for r in rems]  #*
        #for p in ps:
        #    s = [(-n%p + p * (inv[p%w] * (r - -n%p) % w)) / w for r in rems]
            for i, start in enumerate(s):
                slice_size = ((size - start - 1) / p + 1)
                sieve[i + start * w_phi :: p * w_phi] = zero * slice_size
        for p in compress(chain_f(izip(*[count(n+r, w) for r in rems])), sieve):
            yield p

Answer 3

Aquí hay un generador que es un poco más cierto de cómo se hace en Haskell: filtrado contra compuestos de primos conocidos, y luego se agregan los primos restantes a la lista.

def gen_primes():
    primes = []
    i = 2
    while True:
        prime = True
        for p in primes:
            if not (i % p):
                prime = False
                break
        if prime:
            yield i
            primes.append(i)
        i += 1

Answer 4

Aquí hay una implementación complicada basada en heap, que no es mucho más rápida que otras implementaciones basadas en heap (ver la comparación de velocidad en otra respuesta mía), pero usa mucha menos memoria.

Esta implementación usa dos montones (tu y wv), que contienen los mismos elementos numéricos. Cada elemento es un par int. Para encontrar todos los primos hasta q**2 (donde q es un primo), cada montón contendrá como máximo 2*pi(q-1) elementos, donde pi(x) es el número de primos positivos no mayores que x . Entonces el número total de enteros es como máximo 4*pi(floor(sqrt(n))) . (Podríamos ganar un factor de 2 en la memoria empujando la mitad de las cosas al montón, pero eso haría que el algoritmo fuera más lento).

Otros enfoques basados ​​en dict y heap (p. Ej., Erat2b y heap_prime_gen_squares y heapprimegen) almacenan enteros `2 * pi (n) '', porque extienden su montón o dict cada vez que encuentran un primo. A modo de comparación: para encontrar los números primos 1_000_000, esta implementación almacena menos de 4141 enteros, otras implementaciones almacenan más de 1_000_000 enteros.

import heapq

def heap_prime_gen_smallmem():
    yield 2
    yield 3
    f = 5
    fmar3 = 2
    q = 7
    q6 = 7 * 6
    qmar3 = 4
    tu = [(25, 30), (35, 30)]
    vw = [(25, 30), (35, 30)]
    while True:
        qmar3 += 2   
        if qmar3 == 6:  
            qb = q + 4
            q6b = q6 + 24
            qmar3 = 2
        else:
            qb = q + 2
            q6b = q6 + 12
        if q < tu[0][0]:
            d = q * q
            while f < d:
                a, b = vw[0]
                if f < a: 
                    yield f   
                else:
                    a, b = vw[0]
                    heapq.heapreplace(vw, (a + b, b))
                    a, b = vw[0]
                    while f >= a:
                        heapq.heapreplace(vw, (a + b, b))
                        a, b = vw[0]   
                fmar3 += 2
                if fmar3 == 6:
                    f += 4
                    fmar3 = 2
                else:
                    f += 2
            c = q * qb   
            heapq.heappush(tu, (d, q6))
            heapq.heappush(tu, (c, q6))
            heapq.heappush(vw, (d, q6))
            heapq.heappush(vw, (c, q6))
        else:
            a, b = tu[0]
            heapq.heapreplace(tu, (a + b, b))
            a, b = tu[0]  
            while q >= a:
                heapq.heapreplace(tu, (a + b, b))
                a, b = tu[0]
        q = qb
        q6 = q6b

Answer 5

Aquí hay una simple pero no terriblemente lenta usando un montón en lugar de un dict:

import heapq

def heap_prime_gen_squares(): 
    yield 2  
    yield 3  
    h = [(9, 6)]
    n = 5
    while True:
        a, b = h[0]
        while n < a:
            yield n
            heapq.heappush(h, (n * n, n << 1))
            n += 2
        heapq.heapreplace(h, (a + b, b))  # Replace h[0], which is still (a, b).

Mis medidas de velocidad del tiempo de usuario para los primeros 1 millón de primos (los números más pequeños son mejores):

• postponed_sieve (basado en dict): 8.553s
• erat2b (basado en el dictado): 9.513s
• erat2a (basado en el dictado): 10.313s
• heap_prime_gen_smallmem (basado en el montón): 23.935s
• heap_prime_gen_squares (basado en el montón): 27.302s
• heapprimegen (basado en el dictado): 145.029s

Así que los enfoques basados ​​en dict parecen ser los más rápidos.

Answer 6

Dado que el OP solicita una implementación eficiente , aquí hay una mejora significativa del código del estado activo 2002 por David Eppstein / Alex Martelli (visto aquí en su respuesta ): no registre la información de un primo en el diccionario hasta que su casilla se vea entre los candidatos . Reduce la complejidad del espacio por debajo de O (sqrt (n)) en lugar de O (n) , para n primos producidos ( π (sqrt (n log n)) ~ 2 sqrt (n log n) / log (n log n) ~ 2 sqrt (n / log n) ). En consecuencia, la complejidad del tiempo también se mejora, es decir, se ejecuta más rápido .

Crea un "tamiz deslizante" como un diccionario de múltiplos actuales de cada base principal (es decir, por debajo del sqrt del punto de producción actual), junto con sus valores de pasos :

from itertools import count
                                         # ideone.com/aVndFM
def postponed_sieve():                   # postponed sieve, by Will Ness      
    yield 2; yield 3; yield 5; yield 7;  # original code David Eppstein, 
    sieve = {}                           #   Alex Martelli, ActiveState Recipe 2002
    ps = postponed_sieve()               # a separate base Primes Supply:
    p = next(ps) and next(ps)            # (3) a Prime to add to dict
    q = p*p                              # (9) its sQuare 
    for c in count(9,2):                 # the Candidate
        if c in sieve:               # c''s a multiple of some base prime
            s = sieve.pop(c)         #     i.e. a composite ; or
        elif c < q:  
             yield c                 # a prime
             continue              
        else:   # (c==q):            # or the next base prime''s square:
            s=count(q+2*p,2*p)       #    (9+6, by 6 : 15,21,27,33,...)
            p=next(ps)               #    (5)
            q=p*p                    #    (25)
        for m in s:                  # the next multiple 
            if m not in sieve:       # no duplicates
                break
        sieve[m] = s                 # original test entry: ideone.com/WFv4f

(el código original más antiguo aquí fue editado para incorporar los cambios como se ve en la respuesta de Tim Peters , a continuación). ver también this para una discusión relacionada.

El código similar a una rueda 2-3-5-7 funciona ~ 2.15x más rápido (lo cual está muy cerca de la mejora teórica de 3/2 * 5/4 * 7/6 = 2.1875 ).

Answer 7

Este no es originalmente mi código, sin embargo, vale la pena publicarlo. El original se puede encontrar aquí: http://code.activestate.com/recipes/117119/

def gen_primes():
  D = {}
  q = 2  # first integer to test for primality.

  while True:
    if q not in D:
      # not marked composite, must be prime  
      yield q 

      #first multiple of q not already marked
      D[q * q] = [q] 
    else:
      for p in D[q]:
        D.setdefault(p + q, []).append(p)
      # no longer need D[q], free memory
      del D[q]

    q += 1

Es un generador, así que úsalo como cualquier otro.

primes = gen_primes()
for p in primes:
  print p

Se necesitan 1.62 segundos para generar y poner en un conjunto, 1 millón de primos, en mi escritorio.

Answer 8

Esto es similar pero aproximadamente un 5% más rápido que erat2a.

import itertools as it
def erat2b( ):
    D = {  }
    yield 2
    for q in it.islice(it.count(3), 0, None, 2):
        p = D.pop(q, None)
        if p is None:
            D[q*q] = q
            yield q
        else:
            # Only this part is replaced.
            q += p
            while q in D:
                q += p
            D[q] = p

Answer 9

Hace algún tiempo escribí un artículo sobre un generador de primos infinitos:

http://stacktrace.it/2008/01/progetto-eulero-problema-3/

Está en italiano, pero es posible que tengas una traducción molesta con Google: http://tinyurl.com/yzpyeom

Answer 10

Haga un tamiz segmentado , donde el tamaño de un segmento está determinado por la memoria disponible o el tamaño máximo de un conjunto de bits.

Para cada segmento representan los números en algún intervalo [n; n + segment_size) como un conjunto de bits y tamiz con todos los números primos debajo de la raíz cuadrada del límite superior.

Usar un conjunto de bits utiliza menos memoria que una tabla hash o una estructura de datos de árbol, porque está trabajando con conjuntos densos de números.

Answer 11

Otra forma de hacerlo:

import itertools
def primeseq():
    prime = [2]
    num = 0
    yield 2
    for i in itertools.count(3, 2):
        is_prime = True
        for num in prime:
            if i % num == 0:
                is_prime = False
                break
            elif num ** 2 > i: 
                break
        if is_prime:
            prime.append(i)
            yield i

Answer 12

Para la posteridad, aquí hay una reescritura del hermoso algoritmo de Will Ness para Python 3. Se necesitan algunos cambios (los iteradores ya no tienen los métodos .next() , pero hay una nueva función incorporada next() ). Otros cambios son divertidos (usar el nuevo yield from <iterable> reemplaza cuatro declaraciones de yield en el original. Más son para la legibilidad (no soy fanático de usar en exceso ;-) nombres de variable de 1 letra).

Es significativamente más rápido que el original, pero no por razones algorítmicas. La aceleración se debe principalmente a la eliminación de la función add() del original, haciendo eso en línea en su lugar.

def psieve():
    import itertools
    yield from (2, 3, 5, 7)
    D = {}
    ps = psieve()
    next(ps)
    p = next(ps)
    assert p == 3
    psq = p*p
    for i in itertools.count(9, 2):
        if i in D:      # composite
            step = D.pop(i)
        elif i < psq:   # prime
            yield i
            continue
        else:           # composite, = p*p
            assert i == psq
            step = 2*p
            p = next(ps)
            psq = p*p
        i += step
        while i in D:
            i += step
        D[i] = step

Answer 13

Sé que la publicación es antigua, pero me topé con esta pregunta ... El siguiente código se basa en una idea muy simple: un creciente tamiz de Eratóstenes. Esta solución es realmente más lenta que las mejores, pero es fácil de entender y está diseñada para ser legible ...

Usé números enteros para almacenar los resultados del tamiz. En formato binario, un entero es una lista de 0 sy 1 s, 0 en la posición i si i no es un primo, 1 si puede ser un primo. El infinito requerido es el resultado del hecho de que los enteros de Python 3 no tienen límites.

def primes():
    container, size = 1 << 2, 3 # we start with 0b100 (from right to left: 0 and 1 are not primes, 2 is
    last_prime = 1
    while True:
        prime = next((j for j in range(last_prime+1, size) if container & 1 << j), None) # find the next prime
        while not prime:
            container, size = expand(container, size, 2**16) # add 65536 cells and sieve the container
            prime = next((j for j in range(last_prime+1, size) if container & 1 << j), None)
        yield prime
    last_prime = prime

¿Cómo expandir el contenedor? Simplemente agregue un grupo de 1 s a la izquierda del contenedor (en formato binario) y crítelos. Esto es idéntico al tamiz estándar, con una ligera diferencia. En el tamiz estándar, si encontramos un primo i , comenzamos a cruzar las celdas en i*i , con un paso de i .

Aquí, esto puede haberse hecho para la primera parte del contenedor. Solo debemos comenzar al principio de la parte nueva del contenedor si está más lejos que i*i .

def expand(container, size, n):
    new_size = size + n
    container += (1 << (new_size + 1) - 1) - (1 << size) # add n 1''s
    for i in range(2, new_size):
        if container & (1 << i): # i is a prime
            t = sum(1 << j for j in range(max(i, size // i)*i, new_size, i)) # set 1 for all mutiple
            container &= ~t # cross the cells

    return container, new_size

Prueba para un millón de primos:

import itertools
assert 78498 == len(list(itertools.takewhile(lambda p: p<1000000, primes())))

Answer 14

Y otra respuesta, más eficiente en cuanto a la memoria que mi respuesta erat3 aquí:

import heapq

def heapprimegen():
    hp= []
    yield 2
    yield 3
    cn= 3
    nn, inc= 3, 6
    while 1:
        while cn < nn:
            yield cn
            heapq.heappush(hp, (3*cn, 2*cn))
            cn+= 2
        cn= nn+2
        nn, inc= heapq.heappushpop(hp, (nn+inc, inc))

Mantiene un montón (una lista) de múltiplos principales en lugar de un diccionario. Pierde algo de velocidad, obviamente.