Sieve of Erastosthenes highest prime factor in erlang

Question

I have edited the program so that it works(with small numbers) however I do not understand how to implement an accumulator as suggested. The reason why is because P changes throughout the process, therefore I do not know in with which granularity I should break up the mother list. The Sieve of Erastosthenes is only efficient for generating smaller primes, so maybe I should have picked a different algorithm to use. Can anybody recommend a decent algorithm for calculating the highest prime factor of 600851475143? Please do not give me code I would prefer a Wikipedia article of something of that nature.

    -module(sieve).
    -export([find/2,mark/2,primes/1]).

    primes(N) -> [2|lists:reverse(primes(lists:seq(2,N),2,[]))].

    primes(_,bound_reached,[_|T]) -> T;
    primes(L,P,Primes) -> NewList = mark(L,P),
        NewP = find(NewList,P),
        primes(NewList,NewP,[NewP|Primes]).

    find([],_) -> bound_reached;
    find([H|_],P) when H > P -> H;
    find([_|T],P) -> find(T,P). 


    mark(L,P) -> lists:reverse(mark(L,P,2,[])).

    mark([],_,_,NewList) -> NewList;
    mark([_|T],P,Counter,NewList) when Counter rem P =:= 0 -> mark(T,P,Counter+1,[P|NewList]);
    mark([H|T],P,Counter,NewList) -> mark(T,P,Counter+1,[H|NewList]).

I found writing this very difficult and I know there are a few things about it that are not very elegant, such as the way I have 2 hardcoded as a prime number. So I would appreciate any C&C and also advice about how to attack these kinds of problems. I look at other implementations and I have absoulutely no idea how the authors think in this way but its something I would like to master.

I have worked out that I can forget the list up until the most recent prime number found, however I have no idea how I am supposed to produce an end bound (subtle humour). I think there is probably something I can use like lists:seq(P,something) and the Counter would be able to handle that as I use modulo rather than resetting it to 0 each time. Ive only done AS level maths so I have no idea what this is.

I cant even do that can I? because I will have to remove multiples of 2 from the entirety of the list. Im thinking that this algorithm will not work unless I cache data to the harddrive, so I'm back to looking for a better algorithm.

I'm now considering writing an algorithm that just uses a counter and keeps a list of primes which are numbers that do not divide evenly with the previously generated prime numbers is this a good way to do it?

This is my new algorithm that I wrote I think it should work but I get the following error "sieve2.erl:7: call to local/imported function is_prime/2 is illegal in guard" I think this is just an aspect of erlang that I do not understand. However I've no idea how I could find the material to read about it. [Im purposely not using higher order functions etc as I have only read upto the bit on recursion in learnyousomeerlang.org]

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

    primes(N) -> primes(2,N,[2]).

    primes(Counter,Max,Primes) when Counter =:= Max -> Primes;
    primes(Counter,Max,Primes) when is_prime(Counter,Primes) -> primes(Counter+1,Max,[Counter|Primes]);
    primes(Counter,Max,Primes) -> primes(Counter+1,Max,Primes).

    is_prime(X, []) -> true;
    is_prime(X,[H|T]) when X rem H =:= 0 -> false;
    is_prime(X,[H|T]) -> prime(X,T).

The 2nd algorithm does not crash but runs too slowly, I'm thinking that I should reimplement the 1st but this time forget the numbers up until the most recently discovered prime, does anybody know what I could use as an end bound? After looking at other solutions it seems people sometimes just set an arbitrary limit i.e 2 million (this is something I do not really want to do. Others used "lazy" implementations which is what I think I am doing.

Here's someone trying to solve the exact same problem you are, but not in Erlang. The comments and answers may help. http://stackoverflow.com/questions/439814/prime-factor-of-300-000-000-000 — kjw0188, Mar 02 '13 at 09:03

score 4 · Answer 1 · answered Feb 26 '13 at 14:57

This:

lists:seq(2,N div 2)

allocates a list, and as the efficiency guide says, a list requires at least two words of memory per element. (A word is 4 or 8 bytes, depending on whether you have a 32-bit or 64-bit Erlang virtual machine.) So if N is 600851475143, this would require 48 terabytes of memory if I count correctly. (Unlike Haskell, Erlang doesn't do lazy evaluation.)

So you'd need to implement this using an accumulator, similar to what you did with Counter in the mark function. For the stop condition of the recursive function, you wouldn't check for the list being empty, but for the accumulator reaching the max value.

Dmitry Belyaev · Answer 2 · 2013-02-27T06:27:16.727

1

By the way you don't need to test all numbers up to N/2. It is enough to test up to sqrt(N).

Here I wrote a version that takes 20 seconds to find the answer on my machine. It uses kind of lazy list of primes and folding through them. It was fun because I solved some project-euler problems using Haskell quite a long ago and to use the same approach on Erlang was a bit of strange.

edited Feb 27 '13 at 06:27

answered Feb 26 '13 at 18:45

Dmitry Belyaev

2,573
12
17

It would be much faster if you would divide by the factors you find. – starblue Feb 27 '13 at 21:53
Divide what? I have a list of primes, I take each of them up to sqrt(N) and check if a prime is a factor of N. – Dmitry Belyaev Feb 28 '13 at 08:10
1

N, so the remaining number becomes much smaller. – starblue Feb 28 '13 at 17:15
Ah, I see. I find (Prime ^ K) * N1 = N and then can look for factors of N1 instead of N. Thank you! – Dmitry Belyaev Feb 28 '13 at 17:37

score 0 · Answer 3 · answered Feb 28 '13 at 08:14

On your update3:

primes(Counter,Max,Primes) when Counter =:= Max -> Primes;
primes(Counter,Max,Primes) when is_prime(Counter,Primes) -> primes(Counter+1,Max,[Counter|Primes]);
primes(Counter,Max,Primes) -> primes(Counter+1,Max,Primes).

You cannot use your own defined functions as guard clauses as in Haskell. You have to rewrite it to use it in a case statement:

primes(Counter,Max,Primes) when Counter =:= Max ->
    Primes;
primes(Counter,Max,Primes) ->
    case is_prime(Counter,Primes) of
        true ->
            primes(Counter+1,Max,[Counter|Primes]);
        _ ->
            primes(Counter+1,Max,Primes)
    end.

Sieve of Erastosthenes highest prime factor in erlang

3 Answers3