Harsh J

Memoirs of a QWERTY Keyboard

Sieve of Euler in Erlang

9 comments

Learning Erlang is a fun task. Solving PE problems with it even more. You can never really beat PE problems unless you have your own fast implementations of primality tests, prime generation sieves and so on.

A quick way to generate primes, as many would tell you, would be to implement and use the Sieve of Eratosthenes. Its a simple sieve that works by striking off all multiples of primes thus observed in sequence. Wikipedia has an illustration on its working technique.

There is another simple sieve derived from the same idea, known as the Sieve of Euler, which uses a filtering technique based on products of an encountered prime across the remaining elements. This sort of a sieve, again explained in a good manner at Wikipedia, is shorter and easier to implement in Erlang than the traditional Eratosthenes (and I mean proper Eratosthenes, mind you).

The following is the solution I came up with, a couple of days into learning Erlang for great good:

-module(primes).
-export([generatePrimes/1]).

% Sieve of Euler %

generatePrimes (N) when is_integer(N) and (N > 1) ->
    Numbers = lists:seq(2,N),
    BeginIndex = 1,
    generatePrimes(Numbers, [], BeginIndex).

generatePrimes (Numbers, Numbers, _Index) ->
    Numbers;

generatePrimes (Numbers, _PreviousNumbers, Index) ->
    FilteredNumbers = ordsets:subtract(Numbers, multiples(lists:nth(Index, Numbers), Numbers)),
    generatePrimes(FilteredNumbers, Numbers, Index + 1).

multiples (Prime, Numbers) ->
    lists:map(fun(Num) -> Prime * Num end, Numbers).

% Compile as c(primes). %


The function primes:generatePrimes(N) will generate all primes until N. On my 2.93 GHz machine it takes about 2.1s to generate all primes under a million (78498 primes).

As a bonus, you can go ahead and compute single-line inelegant solutions #3 and #10 of Project Euler as:

% Solution to problem 3 %
lists:nth(1, lists:dropwhile(fun(X) -> not(600851475143 rem X == 0) end, lists:reverse(primes:generatePrimes(trunc(math:sqrt(600851475143)))))).

% Solution to problem 10 %
lists:sum(primes:generatePrimes(2000000)).

Notes: There is some redundant multiple calculation involved but I noticed that if I perform a sublist operation to remove that it slows the generation down, so I let it be. The sieving operation also terminates when it notices that no number was filtered in one of its call, as this marks the proper end of it. For instance, to generate all primes under 100, it only sieves through the first five primes.

Update:

Kind commenter Angel (angel at uah dot es) posted a solution (and a multi-process one!) to a proper Eratosthenes’ Sieve (based on the O’Neill’s paper I referenced above), which goes below:

% "Genuine" Erathostenes sieve in erlang
% based on www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf
% 2009 Angel Alvarez
%
% Process  1     2     3     4     5     6     7 ....
%
%         02    03    05    07    11    13    17
% P 02 -> 04
% P 03    04 -> 09
% C 04 == 04    09
% P 05    06    09 -> 25
% C 06 == 06    09    25
% P 07    08    09    25 -> 49
% C 08 == 08    09    25    49
% C 09    10 == 09    25    49
% C 10 == 10    12    25    49
% P 11    12    12    25    49 -> 121
% C 12 == 12    12    25    49    121
% P 13    14    15    25    49    121 -> 169
% C 14 == 14    15    25    49    121    169
% C 15    16 == 15    25    49    121    169
% C 16 == 16    18    25    49    121    169
% P 17    18    18    25    49    121    169 -> 289

-module(multi_erathostenes).
-author("angel at uah dot es").
-compile(export_all).

compare(X,Y) when X > Y -> greater;
compare(X,Y) when X < Y -> smaller;
compare(_,_) -> equal.


worker(Prime,UpperBound)->
	receive
		{test,Number,ParentPID} ->
			ParentPID ! ok,
			case compare(Number,UpperBound) of
				equal ->
%					io:format("[Worker of prime: ~w (~w)] Found composite: ~w~n",[Prime,UpperBound,Number]),
					worker(Prime,UpperBound + Prime);
				smaller ->
					io:format("[Worker of prime: ~w (~w)] Spawn new worker for Prime: ~w~n",[Prime,UpperBound,Number]),
					NextPID =spawn(?MODULE,worker,[Number,Number*Number]),
					worker(Prime,UpperBound,NextPID)
			end;
		stop ->
		    noop
	end.

worker(Prime,UpperBound,NextPID) ->
	receive
		{test,Number,ParentPID} ->
			ParentPID ! ok,
			case compare(Number,UpperBound) of
				equal ->
%					io:format("[Worker of prime: ~w (~w)] Found composite: ~w~n",[Prime,UpperBound,Number]),
					worker(Prime,UpperBound + Prime,NextPID);
				smaller ->
% 					io:format("[Worker of prime: ~w (~w)] passing along: ~w~n",[Prime,UpperBound,Number]),
					NextPID ! {test,Number,self()},
					receive
					    ok -> ok
					end,
					worker(Prime,UpperBound,NextPID);
				greater ->
% 					io:format("[Worker of prime: ~w (~w)] passing along: ~w~n",[Prime,UpperBound,Number]),
					NextPID ! {test,Number,self()},
					receive
					    ok -> ok
					end,
 					worker(Prime,UpperBound+Prime,NextPID)
			end;
		stop ->
		    NextPID ! stop,
		    noop
	end.

start(N) ->
    OnePID=spawn(?MODULE,worker,[2,4]),
    cicle(3,N,OnePID).

cicle(Last,Last,Worker) ->
    Worker ! stop;
cicle(Current,Last,Worker) ->
    Worker ! {test,Current,self()},
    receive
	ok -> ok
    end,
    cicle(Current+1,Last,Worker).

Written by Harsh

June 6th, 2010 at 11:15 pm

9 Responses to 'Sieve of Euler in Erlang'

Subscribe to comments with RSS or TrackBack to 'Sieve of Euler in Erlang'.

  1. Great!!

    If you want to see a JFP-like sieve on erlang i managed to translate it to erlang when i began to learn erlang. I use processes to mimic the lazy features of haskell, thus spawing a new process when needed.

    it uses one erlang process for every sieve cell, in later versions i manage to account and check a variable number of primes on every process, thus reducing message overhead on the first stages of the sieve.

    I hope you like it!

    /angel

    Angel

    8 Jun 10 at 9:26 pm

  2. Its a great implementation, thanks for sharing this!

    P.s. -> I’ve moved your code from comments area into the main article cause it seemed to be breaking the layout. Hope you don’t mind. Ready to move it back if you do or feel that the crediting isn’t visible enough!

    Harsh

    8 Jun 10 at 11:35 pm

  3. Ok the code shows now pretty printed!, No problem for credits.

    You’ll find it a bit slow (message passing is almost as heavier as processes main loop). If you want to optimize this, you should try to accumulate some primes before spawing a new process. Also youll see that messaging is syncronous (see the receive ok -< ok clause) thus limiting a bit the performance on multicore.

    Angel

    9 Jun 10 at 1:20 pm

  4. Hi All:
    Is there a typing error in Angel’s code at line 30 in the “pretty print” and in the source code?

    I think it should read:
    30 compare(X,Y) when X < Y -> smaller;

    Cheers, grandad

    Jim Walls

    14 Jun 10 at 7:43 am

  5. Ah yes, the combination formed a tag and escaped I s’pose. Fixed, thanks!

    Harsh

    14 Jun 10 at 11:39 am

  6. So this is what you were busy with, generating prime numbers :P
    I might try learning erlang (after my currently happening adventures with opengl are over :D ), looks scary though!

    Abhinandh

    15 Jun 10 at 6:44 pm

  7. I’d say go ahead with Haskell instead :)

    Harsh

    16 Jun 10 at 12:04 am

  8. What are your views on Haskell vs Clojure? If a n00b were to be interested in functional programming? Also, recommend a good book.

    Pathik

    23 Jun 10 at 9:06 pm

  9. I’ve not tried/used/compared Clojure, can’t say anything about it. I think Haskell is good. No idea on books, who uses them anymore to learn languages? ;)

    Harsh

    23 Jun 10 at 10:26 pm

Leave a Reply