4
\$\begingroup\$

I just started learning Erlang. Here is my crack at the Seive, based on (https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf):

-module(seive).
-export([seive/1]).

seive(N) ->
    L = lists:seq(1,N),
    {Megass, Ss, Micros} = erlang:timestamp(),
    S = doseive(L, 2),
    {Megase, Se, Microe} = erlang:timestamp(),
    {Megase-Megass, Se-Ss, Microe-Micros, S}.

doseive(L, Index) ->
    Isq = Index*Index-1,
    if Isq > length(L) -> [X || X <- L, X /= -1]; %done crossing out
    true -> 
        case  lists:nth(Index, L) == -1 of
        true -> doseive(L, Index + 1); %nothing to cross this pass
        false -> 
            {L1, L2} = lists:split(Isq, L),
             L3 = lists:map(fun(X) -> if X == -1 orelse (X rem Index == 0 andalso X /= Index) -> -1; true  -> X  end end , L2),
             doseive(lists:append(L1, L3), Index + 1)
    end end.

This takes 52 seconds to generate the first two million primes:

6> seive:seive(2000000).
{0,19,53218,
 [1,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,
  73,79,83,89|...]}

I was expecting it to be faster based on other solutions on the web.

Sticking to serial (I have not even attempted a parallel solution yet, so I'd like to try that on my own first), what are the problems with my code? I looked for "in place" versions of map so that I don't continuously copy lists into new lists, which I think is what is killing me, but the essence of functional is this immutable-ness.


If I was in an imperative language, I would iterate over a hashmap, but there is no "map-update" function for hashmaps. I guess I could try accessing values in a loop, but this solution better lends itself to imperative languages. Would be nice to have something like

mapupdate(F, V)
%updates all values for keys to V where F(key) is true

I also tried a solution with filter:

doseivev2(L, Index) ->
if Index > length(L) -> L;
true ->
    case Index > length(L) of 
        true -> L;
        false ->
            %if you do filter, you dont know what index p^2 is anymore...
            doseivev2(lists:filter(fun(X) -> X =< Index orelse X rem Index /= 0 end, L), Index + 1)
end end.

but it is WAAAAY slower. Does not even finish for 1M.

\$\endgroup\$
1
\$\begingroup\$

You point it, this program copies and traverse the lists too many times. You should think to reduce the size of the list at each prime iteration, you can initialize the first list without all the multiples of 2,3,5 (it is quite easy), you should pass length(L) as a parameter (it needs to traverse the list at each iteration).

Since the sieve algorithm is based on the modification of a set of data, the list is not the best data structure to work with. I think you should have much better result with ets.

[edit] I have tested your version on my laptop using timer:tc/3 and found 24s (your measure is 19.5s, not 52s). a version using ets and the different optimizations I suggest takes 1.2s for the same job.

With 10000000 the times are 252s and 7s, which suggests that with your algorithm the execution time is increasing more than with ets (ration 10 and 6 for a data set 5 time bigger).

Last point, you are returning 1 as a prime which is not correct.

\$\endgroup\$
0
\$\begingroup\$

This pulls together some of the answers and partial answers from the closed question https://stackoverflow.com/questions/146622/sieve-of-eratosthenes-in-erlang/389740, along with performance comparisions. The were all run in the same Ubuntu VM. Hopefully this will be helpful for others who find this question and are looking for a comparison of options.

%% Not really SoE according to https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf
%% Θ(n^2/(log n)^2)
unfaithful_sieve([]) -> [];
unfaithful_sieve([H|T]) -> [H| unfaithful_sieve([ X || X <- T, X rem H /= 0 ])];
unfaithful_sieve(N) -> unfaithful_sieve(lists:seq(2, N)).

3> euler_common:stopwatch(fun() -> euler_common:unfaithful_sieve(2000000) end).
Elapsed time: 325 seconds
ok

Wow, that was slow.

%% SoE according to https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf
%% The main difference is that we clear multiples starting at each Prime^2.
%% Θ(n log log n)
sieve(N) ->
    Array = array:new([{size, N+1}, {fixed,true}, {default, 1}]),
    lists:sublist(
    convert_flags_to_list_of_primes(array:to_list(sieve(Array, N, 2)), 0, []), 3, N
    ).

sieve(Array, Max, Max) -> Array;
sieve(Array, Max, Index) ->
    case array:get(Index, Array) of
    0 -> sieve(Array, Max, Index+1);
    1 -> sieve(clear_multiples(Array, Index), Max, Index+1)
    end.

clear_multiples(Array, Prime) -> clear_multiples(Array, Prime * Prime, Prime).
clear_multiples(Array, Index, Prime) ->
    case Index > array:size(Array) - 1 of
    true -> Array;
    false ->
        clear_multiples(array:set(Index, 0, Array), Index + Prime, Prime)
    end.

convert_flags_to_list_of_primes([], _Index, Primes) -> lists:reverse(Primes);
convert_flags_to_list_of_primes([1|T], Index, Primes) -> convert_flags_to_list_of_primes(T, Index+1, [Index|Primes]);
convert_flags_to_list_of_primes([0|T], Index, Primes) -> convert_flags_to_list_of_primes(T, Index+1, Primes).

4> euler_common:stopwatch(fun() -> euler_common:sieve(2000000) end).
Elapsed time: 4 seconds
ok

I think that one is faithful to SoE, and is obviously much faster. It's awfully long, though.

%% A way better solution from https://stackoverflow.com/questions/146622/sieve-of-eratosthenes-in-erlang
%% (https://stackoverflow.com/users/72346/matt-h)
primes(Prime, Max, Primes, Integers) when Prime > Max -> lists:reverse([Prime|Primes]) ++ Integers;
primes(Prime, Max, Primes, Integers) ->
    [NewPrime|NewIntegers] = [ X || X <- Integers, X rem Prime =/= 0 ],
    primes(NewPrime, Max, [Prime|Primes], NewIntegers).
primes(N) -> primes(2, round(math:sqrt(N)), [], lists:seq(3,N,2)). % skip odds

5> euler_common:stopwatch(fun() -> euler_common:primes(2000000) end).
Elapsed time: 1 seconds
ok

Woah. This is the solution provided by @matt_h in answer to https://stackoverflow.com/questions/146622/sieve-of-eratosthenes-in-erlang/389740#389740. It's faster by a long shot, and it's concise.

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.