Finding Primes Quickly

When it comes to topics like cryptography, prime numbers are absolutely essentially. We often need to determine if a number is prime or generate a list of primes. Let’s look at a few ways this is done.

First of all, recall that we have prime and composite numbers. A prime number, \(p\), is an integer which is only divisible by \(1\) and itself. For example, the numbers \(2\), \(3\), and \(5\) are all prime. There are no numbers that can evenly divide any of these other than \(1\) or itself.

This lends itself to the fact that, other than \(2\), there are no even primes since \(2\) can divide any even number. Any number which is not prime is called composite, because it is composed of products of prime numbers. This is actually defined by something called the Fundamental Theorem of Arithmetic. For example, the number \(9\) is composite, because I can multiply \(3\) by \(3\) and I get \(9\), which also means I can divide \(9\) by \(3\) and get \(3\). This means that any multiples of prime numbers are composite. With \(3\) being prime, \(6\), \(9\), \(12\), and so forth are all composite.

Suppose that we want to determine if a number, \(n\), is prime. We can simply check every number from \(2\) up to \(n-1\) and see if it can evenly divide \(n\). If any number can do this, then \(n\) is not prime. However, this is very inefficient. So, how do we address this?

function isprime(n::Int)::Bool
    n <= 4 && return (n == 2 || n == 3)
    n % 2 == 0 && return false

    for i in 3:2:isqrt(n)
        n % i == 0 && return false
    end

    return true
end

I’ve made a couple small adjustments in the approach above to make things a little more efficient, but we can improve things even more by identifying possible primes up to a point and then checking to see if a number is in that set. Enter the sieve. The idea is that we start from \(2\) and then flip every multiple of \(2\), move to \(3\) and flip every multiple of \(3\), and so forth.

Algorithmically, this is what we want to do:

  1. Create a bit vector where every initial value is \(1\). We’ll manually set \(1\) to \(0\) since \(1\) is not prime.

    Julia indexes from \(1\), but if you are using a language which indexes from \(0\), which is most of them, then set the zeroth component of the bit vector to \(0\) as well.

  2. We iterate over the bit vector and look for the next element equal to \(1\). This would be \(2\) in our example.
  3. Using an inner loop, we’ll iterate through the bit vector and mark every multiple by setting the elements to \(0\) since they are not prime. In our example, this would be \(4\), \(6\), \(8\), \(\dots\).

Once this completes, everything which did not get flipped to \(0\) is a prime. We can now check to see if an element is \(1\) or \(0\) in constant time, \(O(1)\).

function sieve(N::Int)::Vector{Int}
    U = ones(N)
    U[1] = 0

    for p in 2:N
        if U[p] == 1
            for m in p*p:p:N
                U[m] = 0
            end
        end
    end

    U
end

There are probably some additional tweaks we can make to improve efficiency a little bit more, but this should provide a good overview of how we build out a sieve to quickly generate a list of primes for primality checks or other uses.