Haskell Prime Number Generator (part 1 of 3)

Project Euler Problem 7 reads "What is the 10 001st prime number?". Solving this problem requires a prime number generator. I solved this using Haskell, which I'm currently learning, but after I'd done so I learned a lot more by looking at others' solutions. Here I present some of the Haskell prime number generators being discussed in the Project Euler forum, problem 7 thread (I can't link to the forum as you need to have an account and have solved Problem 7 to access that thread):

As I'm just learning Haskell my aim is to annotate these examples at my own level, i.e. for beginners. Note, the examples I'm referring to were meant to solve a particular problem, so in some cases I've tweaked them slightly to make them general prime number generators.

Prime Number Generator from user gauchopuro:

primes (x:xs) = x : primes (filter ((/= 0) . (`mod` x)) xs)

This function must be supplied with a list of the positive integers starting at 2 and it uses the Sieve of Eratosthenes to calculate the list of the primes. For example a function that uses the above to get the first n prime numbers would look like this:

getFirstNPrimes n = take n $ primes [2..]

Or to get the nth prime:

getNthPrime n = primes [2..] !! (n-1)

Although very concise, once you understand the notation involved, primes is also very expressive. Here's what you need to know to parse it:

1. primes uses pattern matching. It expects to be supplied with a list of integers and the(x:xs)to the left of the equals sign means that within the function (i.e. to the right of the equals sign) x will refer to the first element in the list and xs (pronounced exes, as in the plural of x) refers to the remainder of the list.
2. primes is recursive. The list that it generates consists of the first element of the list passed into it and the output of primes when passed the remainder of the list after some filtering.
3. The filter function expects two arguments. The first argument is a predicate and the second is a list. filter returns a new list consisting of the elements from the initial list that yield True when passed to the predicate.
4. The . operator composes two functions i.e. (f . g) x == f (g x).
5. The cons operator ':' appends a new element to the start of a list.

Haskell Lazy Evaluation

Although not strictly required by the primes function shown, it's also useful to understand Haskell's lazy evaluation. It means that Haskell will do just enough work to get to the result it needs. For example, consider the following function:

fib :: Int -> Int -> [Int]
fib a b = a : fib b (a+b)

The first line tells us that this function takes two integers and returns a list of integers. The second line defines the function and says that the list to return is made up of the first integer argument followed by the list produced by calling fib again recursively, with the original second argument as the new first argument, and the sum of the original arguments as the new second argument. This function can be called with 1 1 to generate the Fibonacci sequence as follows:

fib 1 1 = 1 : fib 1 2
        = 1 : 1 : fib 2 3
        = 1 : 1 : 2 : fib 3 5
        = 1 : 1 : 2 : 3 : fib 5 8 etc...

You can see how this generates the Fibonacci sequence, but if you're not familiar with lazy evaluation you'd expect this function to call itself recursively indefinitely and never return. However, lazy evaluation ensures that the function does just enough work to get the result it needs. So if you just want the first 10 Fibonacci numbers you call the function with:

take 10 (fib 1 1)

and the result is:

[1,1,2,3,5,8,13,21,34,55]

After the first 10 elements are generated Haskell has the result it needs and so fib returns without continuing to call recursively.

Understanding primes

Here is primes again:

primes (x:xs) = x : primes (filter ((/= 0) . (`mod` x)) xs)

and, again, it is to be called with the list of integers starting at two, for example like:

primes [2..] !! 101          -- To get the 100th prime

With the above knowledge you can now understand primes as implementing the Sieve of Eratosthenes as follows:

1. It's originally passed a list of integers starting at 2 (2, 3, 4, 5...). The first element of that list, x, will originally be 2. This becomes the first element in the list of prime numbers to be returned.
2. The remainder of the list, xs, is the integers from three upwards. They are filtered with the predicate ((/= 0) . (`mod` x)) to create a new list. Remember x is still 2, so this predicate first finds the remainder when the input is divided by 2, and returns True if that remainder is non-zero. As True is returned for all odd elements in xs, the even numbers are filtered out of the list. This filtered list (3, 5, 7, 9...) is passed into primes recursively.
3. This time the first element in the list is 3, so x  is 3, and this becomes the next prime number that is returned. Now xs is the list of odd integers starting at 5. They are filtered as before with x=3 so that the filtered list is the integers from 5 that don't divide evenly by 2 or 3 (5, 7, 11, 13, 17...). This list is passed recursively into primes.
4. This process continues adding another prime number to the list each time.

This process can be illustrated with the following equivalences:

primes [2, 3, 4, 5, 6...] = 2 : primes [3, 5, 7, 9, 11...]

                          = 2 : 3 : primes [5, 7, 11, 13, 17...]

                          = 2 : 3 : 5 : primes [7, 11, 13, 17, 19...]

                          = 2 : 3 : 5 : 7 : primes [11, 13, 17, 19, 23...]

                          = 2 : 3 : 5 : 7 : 11 : primes [13, 17, 19, 23, 29...]

Although not the most efficient prime number generator, I did find this one to be expressive, and a good learning aid.

I'll follow up in future posts with more prime number generators.

Comparison of C++ and Haskell for Project Euler

I've already solved quite a few of the problems on Project Euler using C++, but as I've recently been learning Haskell, I decided to redo them all using it. Solving real problems is a good way to learn a new language and as I've already done them in C++ I can make some comparisons between the two languages.

Haskell for Maths Problems

Without using Haskell for very long it's easy to see that it is well suited to maths problems. This is because the syntax is optimised for expressing mathematical statements and it often resembles maths notation. For example, consider the following definition for the absolute function, absolute:

absolute n
  | n >= 0    = n
  | otherwise = -n

You don't need to know Haskell to figure out what this function does. Another frequently used example is the factorial function. It can be defined, for example, in the following three ways:

factorial n
  | n == 0    = 1
  | otherwise = n * factorial (n-1)

factorial2 0 = 1
factorial2 n = n * factorial2 (n-1)

factorial3 n = product [1..n]

Not only are these more succinct than in a procedural language, but they also read more naturally.

Haskell for Project Euler, Problem 5

Given Haskell's suitability for maths problems, it should be well suited to Project Euler. I'll demonstrate that here using Problem 5, "What is the smallest positive number that is evenly divisible by all of the numbers from 1 to 20?"

Firstly, here is an explanation of my approach. The solution must be a multiple of all the prime numbers below 20, and as they don't share any common factors the solution must be a multiple of the product of the primes below 20. Therefore I calculate this product and test multiples of it until I find the first one that the non primes below 20 also divide evenly into.

Here is the solution implemented in C++

bool isValidAnswer(int numberToCheck) {
  vector<u64> nonPrimesTo20{ 4, 6, 8, 9, 10, 12, 14, 15, 16, 18 };
  auto isMultiple = true;
  for_each(nonPrimesTo20.begin(), nonPrimesTo20.end(), [&](u64 n) {
    if (numberToCheck%n != 0) {
      isMultiple = false;
    }
  });
  return isMultiple;
}

int main(void) {
  const auto productOfPrimesTo20 = 2*3*5*7*11*13*17*19;
  auto numberToCheck = productOfPrimesTo20;

  while (!isValidAnswer(numberToCheck)) {
    numberToCheck += productOfPrimesTo20;
  }
  cout << "Answer is " << numberToCheck << endl;
  return 0;
}

and here is the same algorithm implemented in Haskell:

productOfPrimesTo20 = 2 * 3 * 5 * 7 * 11 * 13 * 17 * 19
nonPrimesTo20 = [4, 6, 8, 9, 10, 12, 14, 15, 16, 18]
isValidAnswer x = and [mod x n == 0 | n <- nonPrimesTo20]
answer = head [x | x <- [productOfPrimesTo20, productOfPrimesTo20*2..], isValidAnswer x]

Again, it is clear that the Haskell solution is much more succinct than in C++. Not only that, but if you understand the Haskell notation, the Haskell solution is also more expressive, and it's easier to read as you don't need to keep track of a loop or any variables.

If you're not familiar with Haskell then this code will look pretty obfuscated, but I've got posts coming up on Haskell prime number generation that explain many of the techniques used above.