Simple prime number generation in Prolog

prime_gen.pl

/*
  prime_gen (prime_gen.pl) v1.0.0-rc
  Simple prime number generation in Prolog
  Copyright (C) 2023 Julio P. Di Egidio
  http://julio.diegidio.name
  
  This software is released under the GNU-GPLv3 license.
  
  As usual, NO WARRANTY OF ANY KIND is implied.
*/

% Developed and tested on:
% - GNU Prolog 1.5.0 (64 bits)
% - SWI-Prolog (threaded, 64 bits, version 9.0.4)

%% LIBRARY %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%!	prime_gen(+Max:integer, ?P:integer) is nondet.
%	
%	Generates the prime numbers up to `Max`.

prime_gen(Max, _) :-
	\+ integer(Max),
	throw(error(type_error(integer,Max),prime_gen/2)).
prime_gen(Max, P) :-
	prime_gen__do(Max, Ps-Ps, 2, P_), P = P_.

%	prime_gen__do(+Max, +PsD, +N, -P) is nondet.

prime_gen__do(Max, _, N, _) :- Max < N, !, fail.
prime_gen__do(Max, Ps-D, N, P) :-  % ~prime => !, next
	prime_gen__comp(Ps, N), !,
	N1 is N + 1,
	prime_gen__do(Max, Ps-D, N1, P).
prime_gen__do(_, _-[], N, N).      % prime => yield
prime_gen__do(Max, Ps-D, N, P) :-  % prime => store, next
	D = [N|D1],
	N1 is N + 1,
	prime_gen__do(Max, Ps-D1, N1, P).

%	prime_gen__comp(+Ps, +N) is semidet.

prime_gen__comp(Ps, _) :- var(Ps), !, fail.
prime_gen__comp([P|_], N) :-   % P > sqrt(N) => !, fail.
	P * P > N, !, fail.
prime_gen__comp([P|_], N) :-   % found => !, yield.
	0 is N mod P, !.
prime_gen__comp([_|Ps], N) :-  % ~found => next
	prime_gen__comp(Ps, N).

%% WORKBENCH %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%!	test_prime_gen(+Max:integer) is det.
%
%	This is equivalent to:
%	  `test_prime_gen(Max, [print])`
%
%	Example:
%	  `_Max = 100, test_prime_gen(_Max).`

test_prime_gen(Max) :-
	test_prime_gen(Max, [print]).

%!	test_prime_gen(+Max:integer, +Opts:list(callable)) is det.
%
%	Options can be:
%	- `print` : prints the list of generated primes
%	- `rounds(N:integer)` : if `print` not in `Opts`,
%	    averages the elapsed time over `N` tries
%
%	Example:
%	  `_Max = 10000, test_prime_gen(_Max, [rounds(20)]).`

test_prime_gen(Max, Opts) :-
	test_prime_gen__opts(Opts, [Print,Rounds]),
	findall(T,
	(	between(1, Rounds, _),
		test_prime_gen__do(Max, T, Print)
	), Ts),
	(	Rounds == 1
	->	format('~nElapsed CPU time: ~d ms~n', Ts)
	;	test_prime_gen__avgt(Ts, Tm),
		format('~nElapsed CPU time: ~0f ms ', [Tm]),
		format('(avg. over ~d rounds)~n', [Rounds])
	).

test_prime_gen__do(Max, T, Print) :-
	p_cpu_time_init,
	(	Print == true
	->	format('~n', []),
		forall(prime_gen(Max, P), format('~d ', [P])),
		format('~n', [])
	;	forall(prime_gen(Max, _), true)
	),
	p_cpu_time_elapsed(T).

test_prime_gen__avgt(Ts, Tm) :-
	sum_list(Ts, TSum),  % TODO: COMPAT
	length(Ts, TLen),
	Tm is TSum / TLen.

test_prime_gen__opts(Opts, [Print,Rounds]) :-
	(	memberchk(print, Opts)
	->	Print = true
	;	Print = false
	),
	(	Print == false,
		memberchk(rounds(N), Opts),
		integer(N), N > 0
	->	Rounds = N
	;	Rounds = 1
	).

%% COMPATIBILITY %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

:- if(catch(current_prolog_flag(dialect, swi), _, fail)).

p_cpu_time_init :-
	statistics(runtime, _),
	statistics(system_time, _).

p_cpu_time_elapsed(T) :-
	statistics(runtime, [_,Tu]),  % user
	statistics(system_time, [_,Ts]),
	T is Tu + Ts.

:- elif(catch(current_prolog_flag(dialect, gprolog), _, fail)).

p_cpu_time_init :-
	statistics(cpu_time, _).

p_cpu_time_elapsed(T) :-
	statistics(cpu_time, [_,T]).  % user + system

:- else.

:- initialization(
	(	catch((current_prolog_flag(dialect, Dialect),
			  Reason = not_supported(Dialect)
		), _, Reason = not_supported),
		throw(error(system_error(Reason),initialization/1))
	)).

:- endif.