Sieve Of Atkin

16 Mar 2017

Introduction

For some reason, everybody is obsessed with primes. Well, not everybody. Rather mathematicians and someone else. And an interesting problem with primes is how to generate them fast. Really fast. One of the proposals to reach this goal, which is pretty recent, is the Sieve of Atkin.

This sieve has in my opinion too much mathematics concepts for someone who doesn’t have a good mathematical knowledge, like me. That’s why I’m going to try to explain it with “easy” terms.

Paper

If you are interested in this algorithm or the concepts that it covers I beg you to take a look at the paper published by Atkin and Bernstein which is much more complete (although much more complex).

Modular Arithmetic

First we are going to take all the invertible numbers in Z₆₀. Z₆₀ represents the natural numbers in module 60. So Z₆₀ is set of the numbers from 0 to 59.

If a number is >= 60 you can convert it to Z₆₀ by doing the mod 60 operation. And the invertible numbers in Z_n are those numbers that have a modular multiplicative inverse. Which means that if you multiply the number in Z_n by another number in Z_n you obtain 1 in module n. And which numbers have inverse? The numbers in Z_n that are coprime with n.
Two numbers are coprime if their greatest common divisor is 1. For example:

5 is the highest number that both 25 and 60 are divisible by.
The only common divisor of 7 and 60 is 1, so they are coprime.

In the case of Z₆₀, the invertible numbers are:
1 7 11 13 17 19 23 29 31 37 41 43 47 49 53 59

Binary Quadratic Form

Now, we are going to generate numbers with the form n = 60k + s where k is just a natural number and s an invertible number in Z₆₀. If we want to generate all the primes from 60 to a limit, k will take values from 1 to ⌊limit/60⌋-1.

For example we will generate (60+1), (60+7) … (60+59), (120+1) …

We know that all the numbers with the previous form are primes candidates. Each one will be prime if it satisfies some rules which are based on binary quadratic forms. The sieve uses 3 algorithms (named 3.1, 3.2, and 3.3) depending on the module of the s used to generate the number.

There are three cases depending on the form of s. The number can be prime if the equation associated with its form has an odd number of solutions:

• s ∈ 1+4Z Which means that s mod 4 = 1. It will be prime if the equation 4x² + y² = n has an odd number of solutions (x,y) where x,y > 0.
• s ∈ 1+6Z Which means that s mod 6 = 1. It will be prime if the equation 3x² + y² = n has an odd number of solutions (x,y) where x,y > 0.
• s ∈ 11+12Z Which means that s mod 12 = 11. It will be prime if the equation 3x² - y² = n has an odd number of solutions (x,y) where x > y > 0.

For example, the number 67 (60 + 7):

s=7 so it has the form 1 + 6Z. The corresponding equation is 3x² + y² = n. The solutions that can be found are: x=1 y=8. So there are an odd number of solutions. Therefore 67 is still a prime candidate.

The hard part of this is how to find the solutions for those equations. The paper exposes efficient methods for this task but I’m not going to explain them because are quite complex and they are not important to understand the philosophy of this algorithm.

Square Free

Finally, if after the previous algorithm the number is still a prime candidate we have to check if the number is square free. If it is, the number is prime!
By the way, square free means that if you factorize the number you won’t find any factor repeated. For example:

12 is not square free because it’s equal to 2·2·3, and we can see that the factor 2 is repeated.

How do we know if a number is square free? Well, actually we don’t know any algorithm that check it in polynomial time, so you can check if the number is divisible by the square of any of the primes that you have already calculated.

Performance

I’m sure that if you are interested in this algorithm is because of its performance. Unfortunately, I can’t swear that this sieve is faster than the Sieve of Eratosthenes. I don’t pretend to analyze deeply the reasons, although if you are interested in knowing the details I recommend you to have a look at this answer in StackOverflow.

Remember that performance is not only about the number of operations that you have to do to obtain the result; is also about cache misses, instruction level parallelism, branch predictions and much more.

Implementation

I have made an implementation of this algorithm focused on making the code illustrative and simple, forgetting about performance related issues that would make it harder to understand.
Implementation in C++
Implementation in Maple

Acknowledgements

I want to thank Dr. G.M. Diaz-Toca her support, especially for understanding the original paper which requires a good mathematical background. I also want to thank her for her implementation of this sieve in Maple.

References

A.O.L. Atkin, D.J. Bernstein, Prime sieves using binary quadratic forms, Math. Comp. 73 (2004), 1023-1030.