:-main(1000).
main(Max) :- primes(Max, Ps), output(Ps, Os), io:outstream(Os).

output([P|Ps], Os) :- true | Os = [write(P), nl|Os1], output(Ps, Os1).
output([]    , Os) :- true | Os = [].

primes(Max, Ps)    :- true | gen(2, Max, Ns), sift(Ns, Ps).

gen(N0, Max, Ns0)  :- N0 =<Max | Ns0=[N0|Ns1],N1:=N0+1, gen(N1, Max, Ns1) .
gen(N0, Max, Ns0)  :- N0 > Max | Ns0=[] .

sift([P|Xs1], Zs0) :- true | Zs0=[P|Zs1], filter(P, Xs1, Ys0), sift(Ys0, Zs1).
sift([],      Zs0) :- true | Zs0=[] .

filter(P, [X|Xs1],Ys0) :- X mod P=:=0 | filter(P, Xs1, Ys0).
filter(P, [X|Xs1],Ys0) :- X mod P=\=0 | Ys0=[X|Ys1], filter(P, Xs1, Ys1) .
filter(P, [],     Ys0) :- true        | Ys0=[].
