Сито Эратосфена - Поиск простых чисел Python
Просто, чтобы уточнить, это не проблема домашней работы:)
Я хотел найти простые числа для математического приложения, которое я создаю, и наткнулся на Сито Эратосфена.
Я написал реализацию в Python. Но это ужасно медленно. Например, если я хочу найти все простые числа менее 2 миллионов. Это займет > 20 минут. (Я остановил его в этот момент). Как я могу ускорить это?
def primes_sieve(limit):
limitn = limit+1
primes = range(2, limitn)
for i in primes:
factors = range(i, limitn, i)
for f in factors[1:]:
if f in primes:
primes.remove(f)
return primes
print primes_sieve(2000)
UPDATE:
Я закончил делать профилирование этого кода и обнаружил, что довольно много времени было потрачено на удаление элемента из списка. Совершенно понятно, что нужно перебирать весь список (наихудший вариант), чтобы найти элемент, а затем удалить его, а затем скорректировать список (возможно, какая-то копия продолжается?). Во всяком случае, я выбрал список для словаря. Моя новая реализация -
def primes_sieve1(limit):
limitn = limit+1
primes = dict()
for i in range(2, limitn): primes[i] = True
for i in primes:
factors = range(i,limitn, i)
for f in factors[1:]:
primes[f] = False
return [i for i in primes if primes[i]==True]
print primes_sieve1(2000000)
Ответы
Ответ 1
Вы не совсем реализуете правильный алгоритм:
В вашем первом примере primes_sieve
не поддерживает список флагов первичности для удаления/сброса (как в алгоритме), но вместо этого постоянно изменяет размер списка целых чисел, что очень дорого: удаление элемента из списка требует сдвига всех последующих предметы вниз на один.
Во втором примере primes_sieve1
поддерживает словарь флагов первичности, который является шагом в правильном направлении, но он перебирает словарь в неопределенном порядке и избыточно отбирает факторы факторов (а не только факторы простых чисел, как в алгоритм). Вы можете исправить это, отсортировав ключи и пропустив не простые числа (что уже делает это на порядок быстрее), но гораздо эффективнее просто использовать список напрямую.
Правильный алгоритм (со списком вместо словаря) выглядит примерно так:
def primes_sieve2(limit):
a = [True] * limit # Initialize the primality list
a[0] = a[1] = False
for (i, isprime) in enumerate(a):
if isprime:
yield i
for n in range(i*i, limit, i): # Mark factors non-prime
a[n] = False
(Обратите внимание, что это также включает алгоритмическую оптимизацию запуска непростой маркировки в простом квадрате (i*i
) вместо его двойного.)
Ответ 2
def eratosthenes(n):
multiples = []
for i in range(2, n+1):
if i not in multiples:
print (i)
for j in range(i*i, n+1, i):
multiples.append(j)
eratosthenes(100)
Ответ 3
Удаление с начала массива (списка) требует перемещения всех элементов после него. Это означает, что удаление каждого элемента из списка таким образом, начиная с фронта, является операцией O (n ^ 2).
Вы можете сделать это гораздо эффективнее с наборами:
def primes_sieve(limit):
limitn = limit+1
not_prime = set()
primes = []
for i in range(2, limitn):
if i in not_prime:
continue
for f in range(i*2, limitn, i):
not_prime.add(f)
primes.append(i)
return primes
print primes_sieve(1000000)
... или, альтернативно, избегайте переупорядочения списка:
def primes_sieve(limit):
limitn = limit+1
not_prime = [False] * limitn
primes = []
for i in range(2, limitn):
if not_prime[i]:
continue
for f in xrange(i*2, limitn, i):
not_prime[f] = True
primes.append(i)
return primes
Ответ 4
Намного быстрее:
import time
def get_primes(n):
m = n+1
#numbers = [True for i in range(m)]
numbers = [True] * m #EDIT: faster
for i in range(2, int(n**0.5 + 1)):
if numbers[i]:
for j in range(i*i, m, i):
numbers[j] = False
primes = []
for i in range(2, m):
if numbers[i]:
primes.append(i)
return primes
start = time.time()
primes = get_primes(10000)
print(time.time() - start)
print(get_primes(100))
Ответ 5
Я понимаю, что это действительно не отвечает на вопрос о том, как быстро генерировать простые числа, но, возможно, некоторые из них найдут эту альтернативу интересной: поскольку python обеспечивает ленивую оценку через генераторы, сито eratosthenes может быть реализовано точно так, как указано:
def intsfrom(n):
while True:
yield n
n += 1
def sieve(ilist):
p = next(ilist)
yield p
for q in sieve(n for n in ilist if n%p != 0):
yield q
try:
for p in sieve(intsfrom(2)):
print p,
print ''
except RuntimeError as e:
print e
Блок try существует, потому что алгоритм работает до тех пор, пока он не ударит стек и без
попробуйте заблокировать обратную трассировку, нажав фактический вывод, который вы хотите видеть на экране.
Ответ 6
Объединив вклады от многих энтузиастов (включая Гленна Мейнарда и MrHIDEn из комментариев выше), я придумал следующий фрагмент кода в python 2:
def simpleSieve(sieveSize):
#creating Sieve.
sieve = [True] * (sieveSize+1)
# 0 and 1 are not considered prime.
sieve[0] = False
sieve[1] = False
for i in xrange(2,int(math.sqrt(sieveSize))+1):
if sieve[i] == False:
continue
for pointer in xrange(i**2, sieveSize+1, i):
sieve[pointer] = False
# Sieve is left with prime numbers == True
primes = []
for i in xrange(sieveSize+1):
if sieve[i] == True:
primes.append(i)
return primes
sieveSize = input()
primes = simpleSieve(sieveSize)
Время, затраченное на вычисление на моей машине для разных входов мощностью 10, равно:
- 3: 0,3 мс
- 4: 2,4 мс
- 5: 23 мс
- 6: 0,26 с
- 7: 3.1 s
- 8: 33 s
Ответ 7
Простая скорость взлома: когда вы определяете переменную "простые числа", установите шаг 2, чтобы автоматически пропустить все четные числа, и установите начальную точку на 1.
Затем вы можете дополнительно оптимизировать вместо я для простых чисел, использовать для я в простых числах [: round (len (primes) ** 0.5)]. Это значительно увеличит производительность. Кроме того, вы можете исключить числа, заканчивающиеся на 5, для дальнейшего увеличения скорости.
Ответ 8
Используя метод filter
для просеивания списка чисел, можно сделать следующее.
from math import sqrt
def eratosthenes(limit):
lst = range(1, limit)
for i in range(2, int(sqrt(limit)) + 1):
lst = filter(lambda x: x == i or x % i, lst) # sieve
return lst
print eratosthenes(2000000)
Ответ 9
Моя реализация:
import math
n = 100
marked = {}
for i in range(2, int(math.sqrt(n))):
if not marked.get(i):
for x in range(i * i, n, i):
marked[x] = True
for i in range(2, n):
if not marked.get(i):
print i
Ответ 10
Здесь версия, немного более эффективная для памяти (и: правильное сито, а не пробные подразделения). В принципе, вместо того, чтобы хранить массив всех чисел и вычеркивать те, которые не являются первичными, в нем хранится массив счетчиков - по одному для каждого найденного штриха - и проскакивает их впереди предполагаемого простого числа. Таким образом, он использует хранилище, пропорциональное количеству простых чисел, а не до самого высокого значения.
import itertools
def primes():
class counter:
def __init__ (this, n): this.n, this.current, this.isVirgin = n, n*n, True
# isVirgin means it never been incremented
def advancePast (this, n): # return true if the counter advanced
if this.current > n:
if this.isVirgin: raise StopIteration # if this is virgin, then so will be all the subsequent counters. Don't need to iterate further.
return False
this.current += this.n # pre: this.current == n; post: this.current > n.
this.isVirgin = False # when it gone, it gone
return True
yield 1
multiples = []
for n in itertools.count(2):
isPrime = True
for p in (m.advancePast(n) for m in multiples):
if p: isPrime = False
if isPrime:
yield n
multiples.append (counter (n))
Вы заметите, что primes()
является генератором, поэтому вы можете сохранить результаты в списке или использовать их напрямую. Здесь первые n
простые числа:
import itertools
for k in itertools.islice (primes(), n):
print (k)
И, для полноты, здесь таймер для измерения производительности:
import time
def timer ():
t, k = time.process_time(), 10
for p in primes():
if p>k:
print (time.process_time()-t, " to ", p, "\n")
k *= 10
if k>100000: return
На всякий случай, когда вам интересно, я также написал primes()
как простой итератор (используя __iter__
и __next__
), и он работал почти с той же скоростью. Удивили меня тоже!
Ответ 11
Я предпочитаю NumPy из-за скорости.
import numpy as np
# Find all prime numbers using Sieve of Eratosthenes
def get_primes1(n):
m = int(np.sqrt(n))
is_prime = np.ones(n, dtype=bool)
is_prime[:2] = False # 0 and 1 are not primes
for i in range(2, m):
if is_prime[i] == False:
continue
is_prime[i*i::i] = False
return np.nonzero(is_prime)[0]
# Find all prime numbers using brute-force.
def isprime(n):
''' Check if integer n is a prime '''
n = abs(int(n)) # n is a positive integer
if n < 2: # 0 and 1 are not primes
return False
if n == 2: # 2 is the only even prime number
return True
if not n & 1: # all other even numbers are not primes
return False
# Range starts with 3 and only needs to go up the square root
# of n for all odd numbers
for x in range(3, int(n**0.5)+1, 2):
if n % x == 0:
return False
return True
# To apply a function to a numpy array, one have to vectorize the function
def get_primes2(n):
vectorized_isprime = np.vectorize(isprime)
a = np.arange(n)
return a[vectorized_isprime(a)]
Проверьте вывод:
n = 100
print(get_primes1(n))
print(get_primes2(n))
[ 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97]
[ 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97]
Сравните скорость сита эратосфенов и грубой силы на ноутбуке Jupyter. Сито Эратосфена в 539 раз быстрее, чем грубая сила для миллионов элементов.
%timeit get_primes1(1000000)
%timeit get_primes2(1000000)
4.79 ms ± 90.3 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
2.58 s ± 31.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
Ответ 12
Я подумал, что можно просто использовать пустой список в качестве завершающего условия для цикла и придумал это:
limit = 100
ints = list(range(2, limit)) # Will end up empty
while len(ints) > 0:
prime = ints[0]
print prime
ints.remove(prime)
i = 2
multiple = prime * i
while multiple <= limit:
if multiple in ints:
ints.remove(multiple)
i += 1
multiple = prime * i
Ответ 13
import math
def sieve(n):
primes = [True]*n
primes[0] = False
primes[1] = False
for i in range(2,int(math.sqrt(n))+1):
j = i*i
while j < n:
primes[j] = False
j = j+i
return [x for x in range(n) if primes[x] == True]
Ответ 14
Моя самая быстрая реализация:
isprime = [True]*N
isprime[0] = isprime[1] = False
for i in range(4, N, 2):
isprime[i] = False
for i in range(3, N, 2):
if isprime[i]:
for j in range(i*i, N, 2*i):
isprime[j] = False
Ответ 15
я думаю, что это самый короткий код для поиска простых чисел с помощью метода Эратосфена
def prime(r):
n = range(2,r)
while len(n)>0:
yield n[0]
n = [x for x in n if x not in range(n[0],r,n[0])]
print(list(prime(r)))
Ответ 16
Используя немного numpy
, я мог найти все простые числа менее 100 миллионов всего за 2 секунды.
Есть две ключевые особенности, которые следует отметить
- Вырезать кратные
i
только для i
до корня n
- Установка кратных значений от
i
до False
с помощью x[2*i::i] = False
намного быстрее, чем явный цикл python for.
Эти два значительно ускоряют ваш код. Для пределов ниже одного миллиона ощутимое время работы отсутствует.
import numpy as np
def primes(n):
x = np.ones((n+1,), dtype=np.bool)
x[0] = False
x[1] = False
for i in range(2, int(n**0.5)+1):
if x[i]:
x[2*i::i] = False
primes = np.where(x == True)[0]
return primes
print(len(primes(100_000_000)))