I wanted to test my performance optimizing skills, and so wanted to find how fast I could get the first \$n\$ prime numbers. I limited myself to only the standard library, as I'm sure numpy
or another library written in C has a prime generator way faster than Python ever will be - and offloading to a library doesn't really improve my ability to improve performance.
I've implemented both the Sieve of Eratosthenes and the Sieve of Sundaram. The SoE was based off my answer here, and the SoS was based off Wikipedias definition. These are available at the end of the answer.
I improved performance by:
SoE: Vectorizing the creation of primes.
primes[base*2:limit:base] = [False]*(ceil(limit / base) - 2)
SoE: Change the start of the slice from \$2b\$ to \$b^2\$. [1]
primes[base*base:limit:base] = [False]*((((limit - base*base) - 1) // base) + 1)
SoE: Simplify the calculations - addition seems to be faster than multiplication.
primes[base * base::base] = [False] * ((limit - 1) // base - base + 1)
SoE: Use
itertools.compress
, rather than a comprehension.SoS: Vectorize the inner loop.
start = 1 + 3*j step = 1 + 3*j primes[start::step] = [False] * ceil((n - start) / step)
SoS: Vectorize the creation of values that only have one value in the sequence.
When \$\frac{n - \text{start}}{\text{stop}} = \frac{n - (1 + 3j)}{1 + 2j} \le 1\$ is equivalent to \$n \le 2 + 5j\$ we know we can stop at \$j = \frac{n - 2}{5}\$.
multi_stop = (n - 2) // 5 for j in range(1, multi_stop): start = 1 + 3*j step = 1 + 2*j primes[start::step] = [False] * ceil((n - start) / step) if multi_stop >= 1: single_start = multi_stop * 3 + 1 primes[single_start::3] = [False] * ceil((n - single_start) / 3)
- SoS: It doesn't seem like you need the
if
created above, and so you can just save wasted cycles.
I tried defining false = [False]*limit
and slice it, but found it to be slower than creating new lists in the loop.
This got the following prime sieves:
from math import ceil
from itertools import compress
def sieve_eratosthenes(limit):
if limit <= 1:
return []
primes = [True] * limit
for base in range(2, int(limit**0.5 + 1)):
if primes[base]:
primes[base * base::base] = [False] * ((limit - 1) // base - base + 1)
primes[0] = primes[1] = False
return list(compress(range(limit), primes))
def sieve_sundaram(limit):
if limit <= 1:
return []
n = (limit - 1) // 2
primes = [True] * n
for j in range(1, (n - 2) // 5):
start = 1 + 3*j
step = 1 + 2*j
primes[start::step] = [False] * ceil((n - start) / step)
return [2] + [2*i + 1 for i, p in enumerate(primes) if p][1:]
Both are faster than the both the original functions.
Code to generate graphs:
from math import ceil
from itertools import compress
import numpy as np
import matplotlib.pyplot as plt
from graphtimer import Plotter, MultiTimer
def sieve_eratosthenes_orig(limit):
if limit <= 1:
return []
primes = [True] * limit
for base in range(2, int(limit**0.5 + 1)):
if primes[base]:
for composite in range(base * 2, limit, base):
primes[composite] = False
return [num for num, is_prime in enumerate(primes) if is_prime][2:]
def sieve_eratosthenes(limit):
if limit <= 1:
return []
primes = [True] * limit
for base in range(2, int(limit**0.5 + 1)):
if primes[base]:
primes[base * base::base] = [False] * ((limit - 1) // base - base + 1)
primes[0] = primes[1] = False
return list(compress(range(limit), primes))
def sieve_sundaram_orig(limit):
if limit <= 1:
return []
n = (limit - 1) // 2
primes = [True] * n
for j in range(1, n):
for i in range(1, j + 1):
value = i + j + 2*i*j
if value < n:
primes[value] = False
return [2] + [2*i + 1 for i, p in enumerate(primes) if p][1:]
def sieve_sundaram(limit):
if limit <= 1:
return []
n = (limit - 1) // 2
primes = [True] * n
for j in range(1, (n - 2) // 5):
start = 1 + 3*j
step = 1 + 2*j
primes[start::step] = [False] * ceil((n - start) / step)
return [2] + [2*i + 1 for i, p in enumerate(primes) if p][1:]
def sieve_test(limit):
if limit <= 1:
return []
n = (limit - 1) // 2
primes = [True] * n
multi_stop = (n - 2) // 5
for j in range(1, multi_stop):
start = 1 + 3*j
step = 1 + 2*j
primes[start::step] = [False] * ceil((n - start) / step)
return [2] + [2*i + 1 for i, p in enumerate(primes) if p][1:]
def test():
for exp in range(6):
limit = 10 ** exp
assert sieve_test(limit) == sieve_eratosthenes(limit)
def main():
fig, axs = plt.subplots()
axs.set_yscale('log')
axs.set_xscale('log')
(
Plotter(MultiTimer([
sieve_eratosthenes_orig,
sieve_eratosthenes,
sieve_sundaram,
sieve_sundaram_orig,
# sieve_test,
]))
.repeat(5, 5, np.logspace(0.35, 2), args_conv=int)
.min()
.plot(axs, x_label='limit')
)
fig.show()
if __name__ == '__main__':
test()
main()
To use the above code snippet you need to install numpy, matplotlib and graphtimer. All should be available via pypi.
Can they be made faster, or is a different sieve faster?
xrange = range
and ensuring no IndexErrors I get the following graph. Which makes it faster than SoS at 1000 and always slower than SoE. \$\endgroup\$