Back to Questions
948
views6
comments
Stack overflow with Sieve of Eratosthenes
f#

I was trying to write a prime number generator using LazyLists and a Sieve of Eratosthenes, to compare to my other LazyList implementation using traditional iteration / mod checking. The code I came up with is the following:

1
2
3
4
5
6
7
8
9

let rec all_integers_starting_from n = 
    LazyList.consf n (fun () -> all_integers_starting_from (n+1))
and scratch_multiples of n nums = 
    LazyList.filter (fun () -> k%n <> 0) nums
and all_primes_after (LazyList.Cons(p, sieved)) = 
    LazyList.consf p (fun () -> all_primes_after (scratch_multiples_of p sieved))
and primes = 
    let sieved_twos = scratch_multiples_of 2 (all_integers_starting_from 3)
    LazyList.consf 2 (fun () -> all_primes_after sieved_twos)

This code works, but it generates a stack overflow eventually, and I'm not sure why. I don't see an obvious case of non-tail recursion. Anyone have any suggestions?

Thanks

.

How long did it run before getting a stack overflow?

Below is the code that I am running right now with a few typos corrected. For some reason the indenting is not correct...

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26

open Microsoft.FSharp.Collections.LazyList

let rec all_integers_starting_from n = 

    LazyList.consf n (fun () -> all_integers_starting_from (n+1))

and scratch_multiples_of n nums = 

    LazyList.filter (fun k -> k%n <> 0) nums

and all_primes_after (LazyList.Cons(p, sieved)) = 

    LazyList.consf p (fun () -> all_primes_after (scratch_multiples_of p sieved))

and primes = 

    let sieved_twos = scratch_multiples_of 2 (all_integers_starting_from 3)

    LazyList.consf 2 (fun () -> all_primes_after sieved_twos)

    

let all_integers_seq n = n |> all_integers_starting_from |> LazyList.to_seq

Seq.iter (fun a -> printfn "%d" a) (all_integers_seq 5)
By on 2/24/2009 8:57 AM ()

It ran for about 20 seconds before getting a stack overflow. The difference between yours and mine though is you're just printing all integers, not all primes. Try using the sequence "primes" instead of "all_integers"

When I run this code:

1

    Seq.iteri (fun i a -> printfn "p(%d) = %d" i a) (LazyList.to_seq primes)

The last value printed is p(13241) = 142711

By on 2/24/2009 9:10 AM ()

Ah, ok! I'm not very familiar with LazyLists or the ``and'' syntax.

It looks like there is a warning 40 about ``all_primes_after'' in the following code being a recursive object rather than a recursive function. I also get an incomplete pattern match for the Cons in all_primes_after definition, which I think is related.

1
2

(fun () -> all_primes_after sieved_twos)
By on 2/24/2009 9:31 AM ()

Both of those end up being fine I think. The warning 40 means it can't determine at compile time if the access patterns to the LazyList will result in something that isn't infinite. If an error happens it would throw an exception about it at runtime. The incomplete pattern match will never happen, the ones that I didn't match are LazyList.Nil, which is the end of a list, but since the way they've been defined produces a list that NEVER ends, LazyList.Nil would never occur anyway. If it did, I'd get a MatchFailureException instead of a stack overflow. :(

By on 2/24/2009 9:34 AM ()

The problem here is subtle (I investigated in the VS debugger). The internals of LazyList.filter must 'force' the first value to see if it matches its predicate. 'force' executes in a try-with block, and exception blocks prevent tail calls. In this example, you end up with a stack that ends with LazyList.filter filtering out 7s forcing one filtering out 5s, forcing one filtering out 3s, forcing one filtering out 2s. So for each prime you've found so far, you end up with that many 'force' call frames.

I can think of ways around this, but I am trying to decide which way is most elegant. In the meantime, I wanted to share my analysis in case someone else has a good solution.

By on 2/24/2009 10:09 AM ()

I actually suspected this might be the culprit, but wanted to wait until you or someone more knowledgeable replied :P

I thought it might be a decent solution to instead of writing the function "sieve_multiples_of p", make it "sieve_multiples_up_to p primes", and have primes be a list of all primes found so far. But a cursory examination seems to suggest this results in much worse performance than checking primes the traditional way, which is to just check primes up to sqrt(n) for divisibility. In essence that's what the "sieve_multiples_up_to p" does, just in a more subtle manner. Unfortunate, because a true sieve should be quite a bit faster than the normal method of prime checking.

By on 2/24/2009 10:17 AM ()

    Topic tags

    IntelliFactory Offices Copyright (c) 2011-2012 IntelliFactory. All rights reserved.
    Home | Products | Consulting | Trainings | Blogs | Jobs | Contact Us | Terms of Use | Privacy Policy | Cookie Policy     		Built with WebSharper