1%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    2/*	Nan.Numerics.Prime
    3	A simple prime number library
    4	Copyright 2016 Julio P. Di Egidio
    5	<mailto:julio@diegidio.name>
    6	<http://julio.diegidio.name/Projects/Nan.Numerics.Prime/>
    7	
    8	This file is part of Nan.Numerics.Prime.
    9	
   10	Nan.Numerics.Prime is free software: you can redistribute it and/or modify
   11	it under the terms of the GNU General Public License as published by
   12	the Free Software Foundation, either version 3 of the License, or
   13	(at your option) any later version.
   14	
   15	Nan.Numerics.Prime is distributed in the hope that it will be useful,
   16	but WITHOUT ANY WARRANTY; without even the implied warranty of
   17	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   18	GNU General Public License for more details.
   19	
   20	You should have received a copy of the GNU General Public License
   21	along with Nan.Numerics.Prime.  If not, see <http://www.gnu.org/licenses/>.
   22*/

   23%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   24
   25% (SWI-Prolog 7.3.25)
   26
   27:- module(prime_lgc, []).
   28
   29:- public
   30	test_/1,		% +N:posint
   31	div_/2,			% +N:posint, -P:prime
   32	div_rev_/2,		% +N:posint, -P:prime
   33	fact_/2,		% +N:posint, -PFs:list(pfact)
   34	gen_/2,			% +Inf:posint, -P:prime
   35	gen_/3,			% +Inf:posint, +Sup:posint, -P:prime
   36	gen_p_/2,		% +L:prime, -P:prime
   37	gen_p_/3,		% +L:prime, +H:prime, -P:prime
   38	gen_rev_/2,		% +Sup:posint, -P:prime
   39	gen_rev_/3,		% +Inf:posint, +Sup:posint, -P:prime
   40	gen_rev_p_/2,	% +H:prime, -P:prime
   41	gen_rev_p_/3,	% +L:prime, +H:prime, -P:prime
   42	next_/2,		% +N:posint, -P:prime
   43	next_p_/2,		% +P0:prime, -P:prime
   44	prev_/2,		% +N:posint, -P:prime
   45	prev_p_/2,		% +P0:prime, -P:prime
   46	right_/2,		% +N:posint, -P:prime
   47	left_/2.		% +N:posint, -P:prime

A simple prime number library :: logic

To allow for maximum performance, module prime_lgc provides unsafe public (not exported) predicates that user code can call directly instead of calling the safe predicates exported by module prime.

For info on the implementation, see library(nan_numerics_prime).

NOTE: Predicates in this module are unsafe, i.e. do not validate input arguments and are not steadfast.

author
- Julio P. Di Egidio
version
- 1.2.5-beta
See also
- library(nan_numerics_prime)
license
- GNU GPLv3
To be done
- Integrate isqrt function from GMP? */
   68:- use_module(nan_numerics_prime_mem).
   69:- use_module(nan_numerics_prime_whl).
   70:- use_module(nan_numerics_prime_prb).
   71:- use_module(library(lists)).
   72
   73%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 test_(+N:posint) is semidet
True if N is a prime number.
   79test_(N) :-
   80	prime_mem:get_(N), !.
   81test_(N) :-
   82	prime_whl:test_(N, Cert),
   83	test__do(Cert, N).
 div_(+N:posint, -P:prime) is semidet
True if N is a composite number with P its smallest prime divisor.
   89div_(N, P) :-
   90	N > 3,
   91	\+ test_(N),
   92	div__do(N, P).
 div_rev_(+N:posint, -P:prime) is semidet
True if N is a composite number with P its greatest prime divisor.
   98div_rev_(N, P) :-
   99	N > 3,
  100	\+ test_(N),
  101	div_rev__do(N, P).
 fact_(+N:posint, -PFs:list(pfact)) is det
PFs is the list of all prime factors of N in ascending order of the prime divisors.

Elements of PFs are of the form P^F with P the prime divisor and F the corresponding power.

If N is equal to 1 or if N is a prime number, PFs is [N^1].

  113fact_(N, PFs) :-
  114	fact__sel_rev(N, [], PFsRev),
  115	reverse(PFsRev, PFs).
  116
  117fact__sel_rev(N, PFs0, PFs) :-
  118	div_(N, P), !,
  119	N1 is N // P,
  120	fact__pfs(P, PFs0, PFs1),
  121	fact__sel_rev(N1, PFs1, PFs).
  122fact__sel_rev(P, PFs0, PFs) :-
  123	fact__pfs(P, PFs0, PFs).
  124
  125fact__pfs(P, [P^F| PFs], [P^F1| PFs]) :- !, F1 is F + 1.
  126fact__pfs(P, PFs, [P^1| PFs]).
 gen_(+Inf:posint, -P:prime) is multi
 gen_(+Inf:posint, +Sup:posint, -P:prime) is nondet
Generates in ascending order all prime numbers P greater than or equal to Inf, and less than or equal to Sup in the variant with arity 3. Fails if the prime to the left of Sup is less than the prime to the right of Inf.
  136gen_(Inf, P) :-
  137	right_(Inf, L),
  138	gen_p_(L, P).
  139
  140gen_(Inf, Sup, P) :-
  141	gen__lh(Inf, Sup, L, H),
  142	gen_p_(L, H, P).
  143
  144gen__lh(Inf, Sup, L, H) :-
  145	Sup >= Inf,
  146	right_(Inf, L),
  147	left_(Sup, H),
  148	H >= L.
 gen_p_(+L:prime, -P:prime) is multi
 gen_p_(+L:prime, +H:prime, -P:prime) is nondet
Generates in ascending order all prime numbers P starting from L, and up to H in the variant with arity 3. Fails if H is less than L.
  156gen_p_(P, P).
  157gen_p_(L, P) :-
  158	next_p_(L, L1),
  159	gen_p_(L1, P).
  160
  161gen_p_(P, H, P) :- P >= H, !.
  162gen_p_(P, _, P).
  163gen_p_(L, H, P) :-
  164	next_p_(L, L1),
  165	gen_p_(L1, H, P).
 gen_rev_(+Sup:prime, -P:prime) is nondet
 gen_rev_(+Inf:posint, +Sup:posint, -P:prime) is nondet
Generates in descending order all prime numbers P less than or equal to Sup, and greater than or equal to Inf in the variant with arity 3. Fails if Sup is equal to 1 or if the prime to the left of Sup is less than the prime to the right of Inf.
  175gen_rev_(Sup, P) :-
  176	left_(Sup, H),
  177	gen_rev_p_(H, P).
  178
  179gen_rev_(Inf, Sup, P) :-
  180	gen__lh(Inf, Sup, L, H),
  181	gen_rev_p_(L, H, P).
 gen_rev_p_(+H:prime, -P:prime) is multi
 gen_rev_p_(+L:prime, +H:prime, -P:prime) is nondet
Generates in descending order all prime numbers P starting from H, and down to L in the variant with arity 3. Fails if H is less than L.
  190gen_rev_p_(2, 2) :- !.
  191gen_rev_p_(P, P).
  192gen_rev_p_(H, P) :-
  193	prev_p_(H, H1),
  194	gen_rev_p_(H1, P).
  195
  196gen_rev_p_(L, P, P) :- P =< L, !.
  197gen_rev_p_(_, P, P).
  198gen_rev_p_(L, H, P) :-
  199	prev_p_(H, H1),
  200	gen_rev_p_(L, H1, P).
 next_(+N:posint, -P:prime) is det
P is the smallest prime number greater than N.
  206next_(N, P) :-
  207	right_(N, P0),
  208	next__sel(N, P0, P).
  209
  210next__sel(P0, P0, P) :- !,
  211	next_p_(P0, P).
  212next__sel(_, P, P).
 next_p_(+P0:prime, -P:prime) is det
P is the smallest prime number greater than P0.
  218next_p_(P0, P) :-
  219	prime_mem:get_(P0, P), !.
  220next_p_(P0, P) :-
  221	N is P0 + 1,
  222	right_(N, P),
  223	prime_mem:add_(P0, P).
 prev_(+N:posint, -P:prime) is semidet
P is the greatest prime number less than N. Fails if N is less than or equal to 2.
  230prev_(N, P) :-
  231	left_(N, P0),
  232	prev__sel(N, P0, P).
  233
  234prev__sel(P0, P0, P) :- !,
  235	prev_p_(P0, P).
  236prev__sel(_, P, P).
 prev_p_(+P0:prime, -P:prime) is semidet
P is the greatest prime number less than P0. Fails if P is equal to 2.
  243prev_p_(P0, P) :-
  244	prime_mem:get_(P, P0), !.
  245prev_p_(P0, P) :-
  246	N is P0 - 1,
  247	left_(N, P),
  248	prime_mem:add_(P, P0).
 right_(+N:posint, -P:prime) is det
P is the smallest prime number greater than or equal to N.
  254right_(N, P) :-
  255	prime_whl:right_(N, W, Cert),
  256	right__sel(Cert, W, P).
  257
  258right__sel(Cert, P, P) :-
  259	test__do(Cert, P), !.
  260right__sel(_, W, P) :-
  261	W1 is W + 2,
  262	right_(W1, P).
 left_(+N:posint, -P:prime) is semidet
P is the greatest prime number less than or equal to N. Fails if N is equal to 1.
  269left_(N, P) :-
  270	prime_whl:left_(N, W, Cert),
  271	left__sel(Cert, W, P).
  272
  273left__sel(Cert, P, P) :-
  274	test__do(Cert, P), !.
  275left__sel(_, W, P) :-
  276	W1 is W - 2,
  277	left_(W1, P).
  278
  279%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  280
  281%	test__do(+Cert:boolean, +N:posint) is semidet.
  282
  283test__do(true, _) :- !.
  284test__do(_, N) :-
  285	prime_prb:test_(N, _).
  286
  287%	div__do(+N:posint, -P:prime) is semidet.
  288
  289div__do(N, P) :-
  290	% Sup is floor(sqrt(N)),		% TODO: Integrate =isqrt= function from GMP?
  291	gen_p_(2, P),
  292	% (	P > Sup, !, fail			% 
  293	(	P * P > N, !, fail			% 
  294	;	N mod P =:= 0
  295	), !.
  296
  297%	div_rev__do(+N:posint, -P:prime) is semidet.
  298
  299div_rev__do(N, P) :-
  300	Sup is N // 2,
  301	left_(Sup, H),
  302	gen_rev_p_(H, P),
  303	N mod P =:= 0, !.
  304
  305%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%