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:

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 <= y:

            value = x_coef * x**2 + y_coef * y**2
            if value > limit:

            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 sieve

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.0026933066114727716
Number of primes found: 10000
    Total Time: 297.5023943142171
    Time per Iteration: 0.0297502394314171
  • \$\begingroup\$ You might find some interesting information on my question.. THere is a sieve of erasthones and atkins in the answers, and they got some wicked fast times. Its in C#, but you might still get something from it. \$\endgroup\$
    – John Davis
    Jun 11, 2014 at 15:38
  • \$\begingroup\$ Wow! u have coded the algorithm ! Great Man thank you but still it is not producing primes i guess.how to produce primes using this code \$\endgroup\$ Jun 12, 2014 at 5:56
  • \$\begingroup\$ i want to compare speed of both ur code and this one codereview.stackexchange.com/questions/6477/… \$\endgroup\$ Jun 12, 2014 at 6:32

1 Answer 1

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.

  • \$\begingroup\$ You were right to assume I'm using Python 3. If this wasn't a 'for-fun' algorithm, I would be more concerned about backwards-compatibility. Thanks for the input! \$\endgroup\$ Jun 12, 2014 at 12:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.