patternpythonMinor
Sieve of Atkins algorithm
Viewed 0 times
atkinsalgorithmsieve
Problem
I got the inspiration to write this algorithm from this question. The OP cited the Sieve of Atkin as one of the fastest ways to generate primes. So I figured I would give it a try:
The algorithm is based off of pseudocode given by the reference to the Sieve of Atkins above. The reference does say:
This pseudocode is written for clarity. Repeated and wasteful calculations mean that it would run slower than the sieve of Eratosthenes. To improve its efficiency, faster methods must be used to find solutions to the three quadratics.
I wanted to make the pseudocode feel as Pythonic as possible. I also tried to optimize it as much as I could without delving too much into the mathematics.
Do you guys see any improvements, performance or style-wise?
Here's some profiling:
```
import timeit
for x in [100, 1000, 10000]:
time = timeit.timeit('sieve_primes({})'.format(x),
setup='from __main__ import sieve_primes', number=10000)
print('Number of primes found: {}'.format(x))
print('\tTotal time: {}\n\tTime per iteration: {}'.format(time, float(time)/10000))
# Results
Number of primes found: 100
Total Time: 2.3724100288364544
Time per Iteration: 0.00023724100288364545
Number of primes found: 1000
Total Time: 26.933066114727716
Time per Iteration: 0.002693306
def sieve_primes(limit=100):
sieve = [False]*limit
loop_range = range(1, int(limit**(1/2)))
for x, y in [(x,y) for x in loop_range for y in loop_range]:
for x_coef, y_coef, mods in [(4, 1, [1, 5]), (3, 1, [7]), (3, -1, [11])]:
if y_coef == -1 and x limit:
continue
if value % 12 in mods:
sieve[value] = not sieve[value]
for value in loop_range[4:]:
if sieve[value]:
for square in [i * value**2 for i in range(limit//value**2+1)
if i * value**2 < limit]:
sieve[square] = False
return sieveThe algorithm is based off of pseudocode given by the reference to the Sieve of Atkins above. The reference does say:
This pseudocode is written for clarity. Repeated and wasteful calculations mean that it would run slower than the sieve of Eratosthenes. To improve its efficiency, faster methods must be used to find solutions to the three quadratics.
I wanted to make the pseudocode feel as Pythonic as possible. I also tried to optimize it as much as I could without delving too much into the mathematics.
Do you guys see any improvements, performance or style-wise?
Here's some profiling:
```
import timeit
for x in [100, 1000, 10000]:
time = timeit.timeit('sieve_primes({})'.format(x),
setup='from __main__ import sieve_primes', number=10000)
print('Number of primes found: {}'.format(x))
print('\tTotal time: {}\n\tTime per iteration: {}'.format(time, float(time)/10000))
# Results
Number of primes found: 100
Total Time: 2.3724100288364544
Time per Iteration: 0.00023724100288364545
Number of primes found: 1000
Total Time: 26.933066114727716
Time per Iteration: 0.002693306
Solution
for x, y in [(x,y) for x in loop_range for y in loop_range]:Here you could avoid the list creation by using a generator (by changing the
[] to ()) or you could use itertools.product. Probably no difference performance-wise, but the latter is less code.for x_coef, y_coef, mods in [(4, 1, [1, 5]), (3, 1, [7]), (3, -1, [11])]:Again, you are creating lists which takes time. This whole thing you iterate over could be a constant nested tuple that's only allocated once at the top (or even outside the function).
for square in [i * value**2 for i in range(limit//value**2+1)
if i * value**2 < limit]:At this point you can probably guess that I'll recommend doing away with the list. A generator should work. If on Python 2,
xrange is probably better than range as well.(Note, I didn't mention
loop_range. That could be changed into xrange invocations if you are on Python 2, but you seem to be on Python 3.)There's also a question mark I have:
int(limit**(1/2)) – that works in Python 3, but not Python 2, because 1/2 == 0. Since your code seemed to work for you, I assume you are on Python 3. I'm not sure if Python 3 people think making it portable with a dot/cast is a good idea or not (i.e. 1./2). There's always the explicit sqrt if you can't decide.Code Snippets
for x, y in [(x,y) for x in loop_range for y in loop_range]:for x_coef, y_coef, mods in [(4, 1, [1, 5]), (3, 1, [7]), (3, -1, [11])]:for square in [i * value**2 for i in range(limit//value**2+1)
if i * value**2 < limit]:Context
StackExchange Code Review Q#53990, answer score: 3
Revisions (0)
No revisions yet.