Finding primes using a primes map with Haskell and Hugs98

15 Dec 2000

      Hi!

As some of you may know, a Haskell program that prints all the primes can be 
as short as the following:

primes = sieve [2.. ] where
         sieve (p:x) = p : sieve [ n | n <- x, n `mod` p > 0 ]

Now, this program roughly corresponds to the following perl program:

###### SNIP SNIP #####
#!/usr/bin/perl

use strict;

my (@primes, $a, $p);
@primes = (2);
MAIN_LOOP: 
for($a = 3; $a < 1000; $a++)
{
    foreach $p (@primes)
    {
        if ($a % $p == 0)
        {
            next MAIN_LOOP;
        }
    }
    push @primes, $a;
}
print join(", ", @primes);
####### SNIP SNIP #####

The program can be more optimized for both speed and code size, but I wanted
to make it as verbose as possible.

The algorithm keeps a list of the primes, and for each new number checks if it
is divisable by any of them and if not it adds it to the list.

There is a different algorithm which keeps a boolean map which tells whether
the number at that position is prime or not. At start it is initialized to all
trues. The algorithm iterates over all the numbers from 2 to the square root
of the desired bound, and if it encounters a prime number it marks all the
numbers p*p, p*p+p, p*p+2*p, p*p+3*p, etc. as not prime. It is generally
considered a better algorithm than the previous one, because it uses less
costier operations (multiplications and additions instead of modulos.)

The perl program that implements that algorithm is this:

#### SNIP SNIP #####
#!/usr/bin/perl

use strict;

sub primes
{
    my $how_much = shift;

    my (@array, $bound, $a, $b, @primes);

    @array = (1) x $how_much;

    $bound = int(sqrt($how_much))+1;

    for($a=2;$a<=$bound;$a++)
    {
        if ($array[$a])
        {
            for($b=$a*$a;$b<$how_much;$b+=$a)
            {
                $array[$b] = 0;
            }
            push @primes, $a;
        }
    }
    for(;$a<$how_much;$a++)
    {
        if ($array[$a])
        {
            push @primes, $a;
        }
    }

    return @primes;
}

print join(", ", primes(1000));
##### SNIP SNIP ######

Now, I tried writing an equivalent Haskell program and the best I could do was
the following:

---- SNIP SNIP -----
module Primes where

import Prelude
import Array

how_much :: Int
how_much = 1000 

initial_primes_map :: Array Int Bool 
initial_primes_map = array (1, how_much) [ (i,True) | i <- [1 .. how_much] ]

mybound :: Int
mybound = ceiling(sqrt(fromInteger(toInteger(how_much))))

next_primes_map :: Int -> Array Int Bool -> Array Int Bool
next_primes_map a primes_map = 
    if (a == mybound) 
    then primes_map 
    else next_primes_map (a+1) (
        if primes_map!a
        then primes_map // [ (i*a, False) | i <- [a .. (prime_bound a)] ]
        else primes_map
        )

prime_bound :: Int -> Int
prime_bound a = (floor(fromInteger(toInteger(how_much))/fromInteger(toInteger(a))))

get_primes_map :: Array Int Bool
get_primes_map = (next_primes_map 2 initial_primes_map)

list_primes :: Array Int Bool -> Int -> [Int]
list_primes primes_map n = 
    if (n > how_much) 
    then [] 
    else 
    (
        if primes_map!n 
        then n:(list_primes primes_map (n+1)) 
        else list_primes primes_map (n+1)
    )

show_primes = show (list_primes get_primes_map 2)
---- SNIP SNIP -----

The problem is that when running it on hugs98 on a Windows98 computer with
64MB of RAM, I cannot seem to scale beyond 30,000 or so, as my boundary. When
entering how_much as 50,000 I get the following message:

ERROR: Garbage collection fails to reclaim sufficient space

In perl I can scale beyond 100,000, and if I modify the code to use a bit
vector (using vec) to much more. So my question is what am I or hugs are doing
wrong and how I can write better code that implements this specific algorithm.
...
From what I saw I used tail recursion, (and hugs98 has proper tail recursion,
right?), and there's only one primes_map present at each iteration (and thus,
at all), so it shouldn't be too problematic. Does it have to do with the way
hugs98 implements and Int to Bool array?
Regards,

    Shlomi Fish

----------------------------------------------------------------------
Shlomi Fish        shlomif@vipe.technion.ac.il 
Home Page:         http://t2.technion.ac.il/~shlomif/
Home E-mail:       shlomif@techie.com

The prefix "God Said" has the extraordinary logical property of 
converting any statement that follows it into a true one.

Shlomi Fish

Joe English

Adrian Hey

Elke Kasimir

Adrian Hey

tags

participants (4)