This code requires the Multi-Precision Integer Arithmetic package.
/* fermat(N, P, Q) is true if P is the largest factor of N less than or */ /* equal to sqrt(N), and Q = N / P, where N, P and Q are strings */ /* containing integer values, and the value of N is odd and greater than 1.*/ /* This is coded using assert/retract in order to reduce heap usage. */ /* e.g. fermat(`8616460799`, P, Q). */ /* fermat(`999846002329`, P, Q). */ fermat(N0, P0, Q0):- string_bigint(N0, N), less_than([1], N), % n > 1 divide_short(N, 2, _, 1), % n is odd square_root(N, S), % s=floor(sqrt(n)) add(S, S, TwoS), add(TwoS, [1], X), % x=2*s+1 multiply(S, S, S2), subtract(N, S2, R), % r=n-s^2 retractall(f(_,_,_,_,_)), assert(f(R,X,[1],[1],[1])), fermat_1, retract(f(_,_,_,P,Q)), !, bigint_string(P, P0), bigint_string(Q, Q0). fermat_1:- f([0],_,_,_,_), retract(f(_,X,Y,_,_)), !, % Factor found subtract(X, Y, Diff), divide_short(Diff, 2, P, _), % p=(x-y)/2 add(X, Y, Sum), subtract(Sum, [2], Sum1), divide_short(Sum1, 2, Q, _), % q=(x+y-2)/2 assert(f([0],X,Y,P,Q)). fermat_1:- % Factor not yet found retract(f(R,X,Y,P,Q)), !, add(X, [2], X1), subtract(X, R, X_R), fermat_2(X_R, Y, Y1, R1), assert(f(R1,X1,Y1,P,Q)), fermat_1. fermat_2(X_R0, Y0, Y, R):- less_than_or_equal(X_R0, Y0), !, add(Y0, [2], Y), subtract(Y0, X_R0, R). fermat_2(X_R0, Y0, Y, R):- subtract(X_R0, Y0, X_R1), add(Y0, [2], Y1), fermat_2(X_R1, Y1, Y, R).