From: notoria <notoria@users.noreply.github.com>
Date: Wed, 29 Apr 2020 00:44:13 +0000 (+0200)
Subject: Implemented predicate mediants/2 with Stern-Brocot tree
X-Git-Tag: v0.8.123~95^2~2
X-Git-Url: https://git.sagredo.dev/?a=commitdiff_plain;h=1d339f74d14d38b8173e43e6983d336c91f8e8d5;p=scryer-prolog.git

Implemented predicate mediants/2 with Stern-Brocot tree
---

diff --git a/src/prolog/lib/arithmetic.pl b/src/prolog/lib/arithmetic.pl
index 406e8473..ccbdc0dd 100644
--- a/src/prolog/lib/arithmetic.pl
+++ b/src/prolog/lib/arithmetic.pl
@@ -1,4 +1,4 @@
-:- module(arithmetic, [msb/2, lsb/2]).
+:- module(arithmetic, [lsb/2, msb/2, mediants/2]).
 
 :- use_module(library(error)).
 
@@ -25,3 +25,28 @@ msb_(X, M, N) :-
     X1 is X >> 1,
     M1 is M + 1,
     msb_(X1, M1, N).
+
+mediants_precision(1.0e-15).
+
+mediants(Real, Fraction) :-
+    builtins:must_be_number(Real, mediants),
+    Fraction = P/Q,
+    mediants_precision(Eps),
+    R is abs(Real),
+    stern_brocot(Eps, R, P1/Q),
+    P is sign(Real) * P1.
+
+% If 0 <= Eps <= 1e-16 then the search is for "infinite" precision.
+stern_brocot(Eps, Real, Fraction) :-
+    Rn is Real - Eps,
+    Rp is Real + Eps,
+    stern_brocot_(Rn, Rp, 0/1, 1/0, Fraction).
+
+stern_brocot_(Rn, Rp, A/B, C/D, Fraction) :-
+    M0 is A + C,
+    M1 is B + D,
+    M is M0 / M1,
+    (   M < Rn -> stern_brocot_(Rn, Rp, M0/M1, C/D, Fraction)
+    ;   M > Rp -> stern_brocot_(Rn, Rp, A/B, M0/M1, Fraction)
+    ;   Fraction = M0/M1
+    ).