gcl-devel
[Top][All Lists]

## RE: [Gcl-devel] RE: GCL bugs (decode-float, rationalize, expt)

 From: Stavros Macrakis Subject: RE: [Gcl-devel] RE: GCL bugs (decode-float, rationalize, expt) Date: Tue, 7 Jan 2003 13:09:58 -0500

```Camm--

Here is code for rationalize which returns the simplest rational.  I
have done some testing, but of course haven't been able to test against
the new integer-decode-float (see comment in code).

-s

(defun float-rationalize (x)
;;; Given x, returns n/d such that float(n/d)=x and d is minimal.
;;; Assuming integer-decode-float works as stated below, automatically
;;; takes into account the precision of the input float.
(cond ((= x 0) 0)
((< x 0) (- (float-rationalize (- x))))
(t (multiple-value-bind
(mantissa exponent sign)
(integer-decode-float x)
;;; We assume that integer-decode-float reflects the internal
;;; representation exactly.  In particular, the number of digits of
;;; precision in mantissa is exactly the number of digits of
;;; precision in the floating-point number, and not reduced to lowest
;;; terms e.g. (integer-decode-float 1.0) => 2^52 -52 1, and NOT 1 0 1.
;;; Similarly for denormalized numbers: (integer-decode-float
;;; least-positive-double-float) should be 1 -1074 1, and not 2^52
;;; -1126 1.  The CL spec is not clear about all this.
(simplest-in-range
;;; A float represents a half-open interval, but using the closed
;;; interval is OK because the midpoint is always simpler than
;;; the endpoints.  Similarly, it never matters that the interval for
;;; exact powers of 2 (except least-positive-float) is assymmetric.
(/ (- (* mantissa 2) 1)
(expt 2 (- 1 exponent)))
(/ (+ (* mantissa 2) 1)
(expt 2 (- 1 exponent))) )))))

;;; From Appendix C4 of the IEEE Scheme standard
;;; Alan Bawden's implementation of Hardy and Wright's algorithm.
;;; Converted to CL by Stavros Macrakis
(defun simplest-in-range (bottom top)
;; 0 < bottom < top
(multiple-value-bind (bottom-integer bottom-fraction)
(truncate bottom)
(if (= bottom-fraction 0)
bottom
(multiple-value-bind (top-integer top-fraction)
(truncate top)
(if (= bottom-integer top-integer)
(+ bottom-integer
(/ 1 (simplest-in-range
(/ 1 top-fraction)
(/ 1 bottom-fraction))))
(+ 1 bottom-integer))))))

```