In this answer we use clpfd, just like this previous answer did.
:- use_module(library(clpfd)).
For easy head-to-head comparison (later on), we call the predicate presented here n_fac/2
:
n_fac(N_expr,F_expr) :-
N #= N_expr, % eval arith expr
F #= F_expr, % eval arith expr
n_facAux(N,F).
Like in this previous answer, n_fac/2
admits the use of arithmetic expressions.
n_facAux(0,1). % 0! = 1
n_facAux(1,1). % 1! = 1
n_facAux(2,2). % 2! = 2
n_facAux(N,F) :-
N #> 2,
F #> N, % redundant constraint
% to help `n_fac(N,N)` terminate
n0_n_fac0_fac(3,N,6,F). % general case starts with "3! = 6"
The helper predicate n_facAux/2
delegates any "real" work to n0_n_fac0_fac/4
:
n0_n_fac0_fac(N ,N,F ,F).
n0_n_fac0_fac(N0,N,F0,F) :-
N0 #< N,
N1 #= N0+1, % count "up", not "down"
F1 #= F0*N1, % calc `1*2*...*N`, not `N*(N-1)*...*2*1`
F1 #=< F, % enforce redundant constraint
n0_n_fac0_fac(N1,N,F1,F).
Let's compare n_fac/2
and n_factorial/2
!
?- n_factorial(47,F).
F = 258623241511168180642964355153611979969197632389120000000000
; false.
?- n_fac(47,F).
F = 258623241511168180642964355153611979969197632389120000000000
; false.
?- n_factorial(N,1).
N = 0
; N = 1
; false.
?- n_fac(N,1).
N = 0
; N = 1
; false.
?- member(F,[3,1_000_000]), ( n_factorial(N,F) ; n_fac(N,F) ).
false. % both predicates agree
OK! Identical, so far... Why not do a little brute-force testing?
?- time((F1 #\= F2,n_factorial(N,F1),n_fac(N,F2))).
% 57,739,784 inferences, 6.415 CPU in 7.112 seconds (90% CPU, 9001245 Lips)
% Execution Aborted
?- time((F1 #\= F2,n_fac(N,F2),n_factorial(N,F1))).
% 52,815,182 inferences, 5.942 CPU in 6.631 seconds (90% CPU, 8888423 Lips)
% Execution Aborted
?- time((N1 #> 1,N2 #> 1,N1 #\= N2,n_fac(N1,F),n_factorial(N2,F))).
% 99,463,654 inferences, 15.767 CPU in 16.575 seconds (95% CPU, 6308401 Lips)
% Execution Aborted
?- time((N1 #> 1,N2 #> 1,N1 #\= N2,n_factorial(N2,F),n_fac(N1,F))).
% 187,621,733 inferences, 17.192 CPU in 18.232 seconds (94% CPU, 10913552 Lips)
% Execution Aborted
No differences for the first few hundred values of N in 2..sup
... Good!
Moving on: How about the following (suggested in a comment to this answer)?
?- n_factorial(N,N), false.
false.
?- n_fac(N,N), false.
false.
Doing fine! Identical termination behaviour... More?
?- N #< 5, n_factorial(N,_), false.
false.
?- N #< 5, n_fac(N,_), false.
false.
?- F in 10..100, n_factorial(_,F), false.
false.
?- F in 10..100, n_fac(_,F), false.
false.
Alright! Still identical termination properties! Let's dig a little deeper! How about the following?
?- F in inf..10, n_factorial(_,F), false.
... % Execution Aborted % does not terminate universally
?- F in inf..10, n_fac(_,F), false.
false. % terminates universally
D'oh! The first query does not terminate, the second does.
What a speedup! :)
Let's do some empirical runtime measurements!
?- member(Exp,[6,7,8,9]), F #= 10^Exp, time(n_factorial(N,F)) ; true.
% 328,700 inferences, 0.043 CPU in 0.043 seconds (100% CPU, 7660054 Lips)
% 1,027,296 inferences, 0.153 CPU in 0.153 seconds (100% CPU, 6735634 Lips)
% 5,759,864 inferences, 1.967 CPU in 1.967 seconds (100% CPU, 2927658 Lips)
% 22,795,694 inferences, 23.911 CPU in 23.908 seconds (100% CPU, 953351 Lips)
true.
?- member(Exp,[6,7,8,9]), F #= 10^Exp, time(n_fac(N,F)) ; true.
% 1,340 inferences, 0.000 CPU in 0.000 seconds ( 99% CPU, 3793262 Lips)
% 1,479 inferences, 0.000 CPU in 0.000 seconds (100% CPU, 6253673 Lips)
% 1,618 inferences, 0.000 CPU in 0.000 seconds (100% CPU, 5129994 Lips)
% 1,757 inferences, 0.000 CPU in 0.000 seconds (100% CPU, 5044792 Lips)
true.
Wow! Some more?
?- member(U,[10,100,1000]), time((N in 1..U,n_factorial(N,_),false)) ; true.
% 34,511 inferences, 0.004 CPU in 0.004 seconds (100% CPU, 9591041 Lips)
% 3,091,271 inferences, 0.322 CPU in 0.322 seconds (100% CPU, 9589264 Lips)
% 305,413,871 inferences, 90.732 CPU in 90.721 seconds (100% CPU, 3366116 Lips)
true.
?- member(U,[10,100,1000]), time((N in 1..U,n_fac(N,_),false)) ; true.
% 3,729 inferences, 0.001 CPU in 0.001 seconds (100% CPU, 2973653 Lips)
% 36,369 inferences, 0.004 CPU in 0.004 seconds (100% CPU, 10309784 Lips)
% 362,471 inferences, 0.036 CPU in 0.036 seconds (100% CPU, 9979610 Lips)
true.
The bottom line?
- The code presented in this answer is as low-level as you should go: Forget
is/2
!
- Redundant constraints can and do pay off.
- The order of arithmetic operations (counting "up" vs "down") can make quite a difference, too.
- If you want to calculate the factorial of some "large"
N
, consider using a different approach.
- Use clpfd!