2

So I'm making a list of prime numbers to help me learn haskell using simple trial division (no fancy stuff until I get better with the language). I'm trying to use the following code:

primes = 2 : [ x | x <- [3..], all (\p -> (mod x p) /= 0) primes]

This is loaded without an error. However:

>take 2 primes
[2ERROR - C stack overflow

I tried the same thing with nested list comprehensions. It doesn't work. I would guess that I'm making too many recursive calls, but this shouldn't be the case if i'm only computing one prime. In my mind the lazy evaluation should make it so that take 2 primes does something along the lines of:

primes = 2 : [ 3 | all (\p -> (mod 3 p) /= 0) [2] ]

Which doesn't require all that much computation - mod 3 2 == True, so all (\p -> (mod 3 p) /= 0) == True, which means take 2 primes == [2, 3], right? I don't understand why this isn't working. Hopefully someone much more versed in the black magic of functional programming can help me...

This is on HUGS, if that makes any difference.

EDIT- I was able to come up with this solution (not pretty):

primes = 2 : [ x | x <- [3..], all (\p -> (mod x p) /= 0) (takeWhile (<= (ceiling (sqrt (fromIntegral x)))) primes)]

EDIT2- The program works fine when interpreted through HUGS or GHCi, but when I try to compile it with GHC, it outputs test: <<loop>>. Anybody know what the problem is?

JasonMArcher
  • 14,195
  • 22
  • 56
  • 52
Robert Mason
  • 3,949
  • 3
  • 31
  • 44
  • 3
    `all (\p -> (mod x p) /= 0) primes` doesn't terminate because `primes` is an infinite sequence. – Tom Crockett Mar 04 '11 at 21:28
  • but the same idea is used in the canonical list comprehension for the fibonacci sequence (if i'm understanding it correctly) `fibonacci = 1 : 1 : [ a + b | (a, b) <- zip fibonacci (tail fibonacci) ]` – Robert Mason Mar 04 '11 at 21:58
  • the problem is not the use of an infinite list, but the use of `all` with an infinite list... verifying that every element of an arbitrary infinite list satisfies some predicate is not possible – Tom Crockett Mar 04 '11 at 22:05
  • On second thought, instead of using `x <- [3..]`, I should probably use `x <- [ (2 * i) + 1 | i <- [1..] ]` – Robert Mason Mar 04 '11 at 23:00
  • The "best" way to do it is [The Genuine Sieve of Eratosthenes](http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf). But the "trial division" approach (also described in the paper) is usually good enough for playing around with Haskell. – Dan Burton Mar 05 '11 at 18:03

3 Answers3

7

Hugs shouldn't do this, but the code is broken anyway so it doesn't matter. Consider:

primes = 2 : [ x | x <- [3..], all (\p -> (mod x p) /= 0) primes]

How do you determine if 3 is prime? well, does mod 3 2 == 0? No. Does mod 3 ??? == 0? OOPS! What is the next element of primes after two? we don't know, we are trying to compute it. You need to add an ordering constraint that adds 3 (or any other x) once all p elem primes less than sqrt x have been tested.

Thomas M. DuBuisson
  • 64,245
  • 7
  • 109
  • 166
  • I believe you should look more closely... 3 _IS_ prime because that's a does _NOT_ equal – Robert Mason Mar 04 '11 at 21:54
  • 4
    His comment is right. That code will attempt to read *ALL* of primes to determine if 3 is prime. But it hasn't yet computed all of primes when it attempts to determine if 3 is prime, so it recurses into itself. It's infinite recursion, which is giving you the error you're seeing. – Carl Mar 04 '11 at 22:08
  • @Robert I think you should look more closely... I wasn't saying 3 isn't prime I was walking through how the given code would compute the 2nd element in the list to show how the code is broken. – Thomas M. DuBuisson Mar 04 '11 at 22:24
  • Thank you for your help. I apologize, I didn't understand what you meant. I come from C++, and (wrongly) assumed that inside a list comprehension the list would be treated as finite, i.e. while you're building the list, if you ask for the list you get the portion of the list that has finished building. – Robert Mason Mar 04 '11 at 22:53
  • Clearly my answer assumed more than it should. I didn't realize there was this alternate notion of list comprehension. – Thomas M. DuBuisson Mar 05 '11 at 01:04
3

The documentation for all says "For the result to be True, the list must be finite" http://hackage.haskell.org/packages/archive/base/latest/doc/html/Prelude.html#v:all

Tim Perry
  • 3,066
  • 1
  • 24
  • 40
0

The previous answers explained why the original comprehension didn't work, but not how to write one that would work.

Here is a list comprehension that recursively, lazily (albeit not efficiently) computes all primes:

let primes = [x | x <- 2:[3,5..], x == 2 || not (contains (\p -> 0 == (mod x p)) (takeWhile (\b -> (b * b) < x) primes))]

Obviously we don't need to check mod x p for all primes, we only need to do it for primes less than the sqrt of the potential prime. That's what the takeWhile is for. Forgive the (\b -> (b * b) < x) this should be equivalent to (< sqrt x) but the Haskell type system didn't like that.

The x == 2 prevents the takeWhile from executing at all before we've added any elements to the list.

Paul Wheeler
  • 18,988
  • 3
  • 28
  • 41
  • Apparently instead of (b * b) < x you could use (fromIntegral b) < (sqrt (fromIntegral x))) but that looks gross. – Paul Wheeler Jul 11 '11 at 02:25
  • For even more fun with with infinite series and list comprehensions: let twinPrimes = [x | x <- zip primes (tail primes), (fst x) + 2 == (snd x)] – Paul Wheeler Jul 11 '11 at 02:58
  • Simplify simplify: let twinPrimes = [(x, x') | (x, x') <- zip primes (tail primes), x + 2 == x']; – Paul Wheeler Jul 11 '11 at 20:18