Как реализовать эффективный бесконечный генератор простых чисел в Python?
Это не домашнее задание, мне просто любопытно.
Ключевое слово здесь - БЕСКОНЕЧНОСТЬ.
Я хочу использовать это как for p in primes()
. Я считаю, что это встроенная функция в Haskell.
Итак, ответ не может быть таким наивным, как "Просто сделай сито".
Прежде всего, вы не знаете, сколько последовательных простых чисел будет потреблено. Ну, предположим, вы могли бы придумать 100 из них за раз. Вы бы использовали тот же подход сита, а также формулу частоты простых чисел?
Я предпочитаю непараллельный подход.
Спасибо за чтение (и написание ;))!
Переведено автоматически
Ответ 1
“Если бы я видел дальше ...”
erat2
Функция из кулинарной книги может быть дополнительно ускорена (примерно на 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
not (x&1)
Проверка проверяет, что x
является нечетным. Однако, поскольку оба q
и p
являются нечетными, при добавлении 2*p
можно избежать половины шагов вместе с проверкой на нечетность.
erat3
Если не возражать против небольшой дополнительной навороченности, erat2
можно ускорить на 35-40% с помощью следующих изменений (ПРИМЕЧАНИЕ: требуется Python 2.7 + или Python 3 + из-за 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
erat3
Функция использует тот факт, что все простые числа (за исключением 2, 3, 5) по модулю 30 приводят только к восьми числам: тем, которые включены в MODULOS
frozenset . Таким образом, после получения начальных трех простых чисел мы начинаем с 7 и работаем только с кандидатами.
Для фильтрации кандидатов используется itertools.compress
функция; “волшебство” заключается в MASK
последовательности; MASK
имеет 15 элементов (в каждых 30 числах 15 нечетных чисел, как выбрано itertools.islice
функцией) с 1
для каждого возможного кандидата, начиная с 7. Цикл повторяется, как указано в itertools.cycle
функции.
Введение фильтрации кандидатов требует еще одной модификации: or (x%30) not in MODULOS
проверки. erat2
Алгоритм обработал все нечетные числа; теперь, когда erat3
алгоритм обрабатывает только r30 кандидатов, нам нужно убедиться, что все D.keys()
могут быть только такими -ложными— кандидатами.
Тесты
Результаты
На сервере Atom 330 Ubuntu 9.10 версий 2.6.4 и 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
На домашнем сервере AMD Geode LX Gentoo, Python 2.6.5 и 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
Тестовый код
primegen.py
Модуль содержит функции erat2
, erat2a
и erat3
. Далее следует сценарий тестирования:
#!/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
Ответ 2
Поскольку OP запрашивает эффективную реализацию, вот значительное улучшение кода active state 2002 Дэвида Эппштейна / Алекса Мартелли (см. Здесь в его ответе): не записывайте информацию о простом числе в словарь, пока его квадрат не будет виден среди кандидатов. Снижает сложность пространства до уровня ниже O (sqrt(n)) вместо O (n) для получения n простых чисел ( π (sqrt(n log n)) ~ 2 sqrt (n log n) / log(n log n) ~ 2 sqrt(n / log n) ). Следовательно, также повышается временная сложность, т.е. Он выполняется быстрее.
Создает "скользящее сито" в виде словаря текущих кратных каждому базовому простому числу (т.е. Ниже sqrt текущей производственной точки) вместе с их значениями шага:
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
(старый исходный код здесь был отредактирован, чтобы включить изменения, как показано в ответе Тима Питерса, ниже). смотрите также Это для соответствующего обсуждения.
Аналогичный код на колесе 2-3-5-7, основанный на коде , выполняется в ~ 2,15 раза быстрее (что очень близко к теоретическому улучшению 3/2 * 5/4 * 7/6 = 2.1875
).
Обновление 2022 года: Недавно я случайно наткнулся на эту старую штуку "NESL" 1990-х годов, которая на самом деле sqrt
трюк с рекурсией, что и на Python.. Итак, ничто не ново под солнцем. :)
Код может быть напрямую дополнен, чтобы начать перечисление простых чисел непосредственно с заданного значения. Это можно увидеть в этой записи на основе JS.
Ответ 3
Для потомков, вот переписанный прекрасный алгоритм Уилла Несса для Python 3. Необходимы некоторые изменения (у итераторов больше нет .next()
методов, но есть новая next()
встроенная функция). Другие изменения предназначены для развлечения (использование нового yield from <iterable>
заменяет четыре yield
оператора в оригинале. Дополнительные сведения приведены для удобства чтения (я не фанат чрезмерного использования ;-) имен переменных из 1 буквы).
Это значительно быстрее, чем оригинал, но не по алгоритмическим причинам. Ускорение в основном связано с удалением функции add()
оригинала, вместо этого она встроена.
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
Ответ 4
Изначально это не мой код, однако его стоит опубликовать. Оригинал можно найти здесь: 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
Это генератор, поэтому используйте его, как любой другой.
primes = gen_primes()
for p in primes:
print p
Для генерации и помещения в набор из 1 миллиона простых чисел на моем рабочем столе требуется 1,62 секунды.