Sieve Of Atkin Implementation in Python

Question

I am trying to implement the algorithm of Sieve of Atkin given in Wikipedia Link as below:

What I've tried so far is the implementation in Python given by following Code:

import math
is_prime = list()
limit = 100
for i in range(5,limit):
    is_prime.append(False)

for x in range(1,int(math.sqrt(limit))+1):
    for y in range(1,int(math.sqrt(limit))+1):
        n = 4*x**2 + y**2

        if n<=limit and (n%12==1 or n%12==5):
            # print "1st if"
            is_prime[n] = not is_prime[n]
        n = 3*x**2+y**2
        if n<= limit and n%12==7:
            # print "Second if"
            is_prime[n] = not is_prime[n]
        n = 3*x**2 - y**2
        if x>y and n<=limit and n%12==11:
            # print "third if"
            is_prime[n] = not is_prime[n]

for n in range(5,int(math.sqrt(limit))):
    if is_prime[n]:
        for k in range(n**2,limit+1,n**2):
            is_prime[k] = False
print 2,3
for n in range(5,limit):
    if is_prime[n]: print n

Now I get error as

is_prime[n] = not is_prime[n]
IndexError: list index out of range

this means that I am accessing the value in list where the index is greater than length of List. Consider the Condition when x,y = 100, then of-course the condition n=4x^2+y^2 will give value which is greater than length of list. Am I doing something wrong here? Please help.

EDIT 1 As suggested by Gabe, using

is_prime = [False] * (limit + 1)

insted of :

for i in range(5,limit):
    is_prime.append(False)

did solved the problem.

The Code is working fine when added the line as suggested by gabe. — Pant, Feb 15 '14 at 02:53

score 6 · Answer 1 · answered Nov 30 '14 at 17:03

Here is a solution

import math

def sieveOfAtkin(limit):
    P = [2,3]
    sieve=[False]*(limit+1)
    for x in range(1,int(math.sqrt(limit))+1):
        for y in range(1,int(math.sqrt(limit))+1):
            n = 4*x**2 + y**2
            if n<=limit and (n%12==1 or n%12==5) : sieve[n] = not sieve[n]
            n = 3*x**2+y**2
            if n<= limit and n%12==7 : sieve[n] = not sieve[n]
            n = 3*x**2 - y**2
            if x>y and n<=limit and n%12==11 : sieve[n] = not sieve[n]
    for x in range(5,int(math.sqrt(limit))):
        if sieve[x]:
            for y in range(x**2,limit+1,x**2):
                sieve[y] = False
    for p in range(5,limit):
        if sieve[p] : P.append(p)
    return P

print sieveOfAtkin(100)

score 3 · Accepted Answer · edited Feb 17 '14 at 05:03

3

You problem is that your limit is 100, but your is_prime list only has limit-5 elements in it due to being initialized with range(5, limit).

Since this code assumes it can access up to limit index, you need to have limit+1 elements in it: is_prime = [False] * (limit + 1)

Note that it doesn't matter that 4x^2+y^2 is greater than limit because it always checks n <= limit.

edited Feb 17 '14 at 05:03

Jim Ferrans

29,992
12
55
83

answered Feb 14 '14 at 15:32

Gabe

82,547
12
135
231

I added range(5, limit+1) while appending in is_prime and yet receive the error. – Pant Feb 14 '14 at 15:36
1

`range(x, y)` returns `y-x` elements, so your `is_prime` list will be indexed from 0 to 95. You need it to be indexed from 0 to 101. Use `is_prime = [0] * (limit + 1)` to create the list. – Gabe Feb 14 '14 at 15:51
Gabe can you refer.to the link given in wikipedia. It has given range 1 to 95 including both – Pant Feb 14 '14 at 15:55
They initialize their list with `False`, not `0`, so it should probably be `[False] * (limit + 1)`. (Not that it makes a difference with the implementation, but it should make it more similar to the original) – Xymostech Feb 14 '14 at 15:57
@SarvagyaPant: Do you understand that in Python lists always start out at index 0? No matter what number you put in your `range()` call, the very first time you call `is_prime.append()`, you are creating an index 0. The 95th time you call `is_prime.append()`, you are making the list 95 elements long regardless of the fact that `i` is 100! – Gabe Feb 14 '14 at 16:09
So silly of me. I got it. I will check the program again with correction. – Pant Feb 14 '14 at 16:22

berdzi · Answer 3 · 2022-03-17T10:54:48.590

This is the optimized implementation proposed by Zsolt KOVACS:

    import math
    import sys
    
    def sieveOfAtkin(limit):
        P = [2,3]
        r = range(1,int(math.sqrt(limit))+1)
        sieve=[False]*(limit+1)
        for x in r:
            for y in r:
                xx=x*x
                yy=y*y
                xx3 = 3*xx
                n = 4*xx + yy
                if n<=limit and (n%12==1 or n%12==5) : sieve[n] = not sieve[n]
                n = xx3 + yy
                if n<=limit and n%12==7 : sieve[n] = not sieve[n]
                n = xx3 - yy
                if x>y and n<=limit and n%12==11 : sieve[n] = not sieve[n]
        for x in range(5,int(math.sqrt(limit))):
            if sieve[x]:
                xx=x*x
                for y in range(xx,limit+1,xx):
                    sieve[y] = False
        for p in range(5,limit):
            if sieve[p] : P.append(p)
        return P
    
    primes = sieveOfAtkin(int(sys.argv[1]))    
    print (primes)

You pass upper limit as the first argument. This program runs in about 6s on my machine comparing to the original which runs in 21s for 10 million limit. What I did:

replaced exponentiation with multiplication
precalculated some multiplications

Sieve Of Atkin Implementation in Python

3 Answers3

Linked