Project Euler

Prime Number Generator

In preparation for later Project Euler challenges, I have been playing around with a prime number generator. My initial method was to test each test_value, n, against all preceding primes: e.g.:

let n = 20
if \exists x \in [2, 3, 5, 7, 11, 13, 17, 19] , that divides n evenly, then n is not prime

This is a vast improvement on the naive approach, which is to test n against all numbers from 2 through n – 1. However, we can do better by testing only those numbers up to (n / 2), or even better by testing only those numbers up to \sqrt{n}.

Well, thanks to an answer here, I have found we can do even better again by testing only 2, 3, and all values that are a multiple of 6 plus or minus one. This uses the observation that all primes can be expressed as (6k +1) or (6k – 1). This is a huge improvement! So here is a prime number generator in python:

import math

def is_prime(n):
    return all((n % i) for i in six_k(n))

def six_k(n):
    ks = [2, 3]
    ks.extend(i for i in range(5, int(math.sqrt(n) + 1))\
              if not (i + 1) % 6 or not (i - 1) % 6)
    return ks

def first_k_primes(n, k, my_primes):
    if len(my_primes) >= k:
        return my_primes

    if is_prime(n):
        my_primes.append(n)

    return first_k_primes(n + 2, k, my_primes)

There are 3 parts to this generator: first_k_primes() generates a list of the first k number of primes, e.g. the first 5 primes are [2, 3, 5, 7, 11]. It calls the is_prime() function which returns whether a test-case (n) is prime or not, which it does by testing (n) against all the primes up to the \sqrt{n}.

In the six_k() function, primes are generated on the fly, by testing against the modulo of (6k +- 1). It is expensive to calculate the primes in this way, on the fly – but I believe the combination of the (6k +- 1) and the generator expression inside ks.extend() help to run this in the best possible time. In fact,  I ran this script with the cProfile utility for the first 200 primes, and six_k() took a full 50 % of the total running time (though it still completed in 0.013 seconds). The alternative is to memoize the list of primes – but this takes up memory and is probably best for very large values of (n), or when you are running this script over and over again. I am hoping, that for most Project Euler challenges, this will be much faster than brute force and therefore fast enough.

(should also note that recursively calculating the first_k_primes() is not a particularly good option in Python! The native recursion depth means it maxs out at approx. the 200th prime number. I will have to tweak this to use a generator expression instead.)

I don’t believe this is the simplest prime number generator around, the six_k() function adds an extra level of complexity. But this is the price of optimisation. If you could always optimise your running time, while also reducing the LOC, then we would probably write optimal algorithms to being with!

This prime number generator will be saved in a library file (myprimes.py), expect to see it referenced in future Project Euler challenges 🙂

— EDIT —
Revised first_k_primes() is below, with a using a generator to yield successive primes:

def first_k_primes(k):
    k, n = k - 2, 5
    while k:
        if is_prime(n):
            yield n
            k = k -1
        n = n + 2

if __name__ == '__main__':
    my_primes = [2, 3]
    my_primes.extend(i for i in first_k_primes(200))
    print(my_primes)

No longer constrained by the recursion limit.

Advertisements

3 thoughts on “Prime Number Generator

  1. Pingback: Haskell Prime Number Generator | :: NickBurns

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s