1 ;;;; This file contains the definitions of float-specific number
2 ;;;; support (other than irrational stuff, which is in irrat.) There is
3 ;;;; code in here that assumes there are only two float formats: IEEE
4 ;;;; single and double. (LONG-FLOAT support has been added, but bugs
5 ;;;; may still remain due to old code which assumes this dichotomy.)
7 ;;;; This software is part of the SBCL system. See the README file for
10 ;;;; This software is derived from the CMU CL system, which was
11 ;;;; written at Carnegie Mellon University and released into the
12 ;;;; public domain. The software is in the public domain and is
13 ;;;; provided with absolutely no warranty. See the COPYING and CREDITS
14 ;;;; files for more information.
16 (in-package "SB!KERNEL")
20 (eval-when (:compile-toplevel :load-toplevel :execute)
22 ;;; These functions let us create floats from bits with the
23 ;;; significand uniformly represented as an integer. This is less
24 ;;; efficient for double floats, but is more convenient when making
25 ;;; special values, etc.
26 (defun single-from-bits (sign exp sig)
27 (declare (type bit sign) (type (unsigned-byte 24) sig)
28 (type (unsigned-byte 8) exp))
30 (dpb exp sb!vm:single-float-exponent-byte
31 (dpb sig sb!vm:single-float-significand-byte
32 (if (zerop sign) 0 -1)))))
33 (defun double-from-bits (sign exp sig)
34 (declare (type bit sign) (type (unsigned-byte 53) sig)
35 (type (unsigned-byte 11) exp))
36 (make-double-float (dpb exp sb!vm:double-float-exponent-byte
38 sb!vm:double-float-significand-byte
39 (if (zerop sign) 0 -1)))
40 (ldb (byte 32 0) sig)))
41 #!+(and long-float x86)
42 (defun long-from-bits (sign exp sig)
43 (declare (type bit sign) (type (unsigned-byte 64) sig)
44 (type (unsigned-byte 15) exp))
45 (make-long-float (logior (ash sign 15) exp)
46 (ldb (byte 32 32) sig)
47 (ldb (byte 32 0) sig)))
53 (defconstant least-positive-single-float (single-from-bits 0 0 1))
54 (defconstant least-positive-short-float least-positive-single-float)
55 (defconstant least-negative-single-float (single-from-bits 1 0 1))
56 (defconstant least-negative-short-float least-negative-single-float)
57 (defconstant least-positive-double-float (double-from-bits 0 0 1))
59 (defconstant least-positive-long-float least-positive-double-float)
60 #!+(and long-float x86)
61 (defconstant least-positive-long-float (long-from-bits 0 0 1))
62 (defconstant least-negative-double-float (double-from-bits 1 0 1))
64 (defconstant least-negative-long-float least-negative-double-float)
65 #!+(and long-float x86)
66 (defconstant least-negative-long-float (long-from-bits 1 0 1))
68 (defconstant least-positive-normalized-single-float
69 (single-from-bits 0 sb!vm:single-float-normal-exponent-min 0))
70 (defconstant least-positive-normalized-short-float
71 least-positive-normalized-single-float)
72 (defconstant least-negative-normalized-single-float
73 (single-from-bits 1 sb!vm:single-float-normal-exponent-min 0))
74 (defconstant least-negative-normalized-short-float
75 least-negative-normalized-single-float)
76 (defconstant least-positive-normalized-double-float
77 (double-from-bits 0 sb!vm:double-float-normal-exponent-min 0))
79 (defconstant least-positive-normalized-long-float
80 least-positive-normalized-double-float)
81 #!+(and long-float x86)
82 (defconstant least-positive-normalized-long-float
83 (long-from-bits 0 sb!vm:long-float-normal-exponent-min
84 (ash sb!vm:long-float-hidden-bit 32)))
85 (defconstant least-negative-normalized-double-float
86 (double-from-bits 1 sb!vm:double-float-normal-exponent-min 0))
88 (defconstant least-negative-normalized-long-float
89 least-negative-normalized-double-float)
90 #!+(and long-float x86)
91 (defconstant least-negative-normalized-long-float
92 (long-from-bits 1 sb!vm:long-float-normal-exponent-min
93 (ash sb!vm:long-float-hidden-bit 32)))
95 (defconstant most-positive-single-float
96 (single-from-bits 0 sb!vm:single-float-normal-exponent-max
97 (ldb sb!vm:single-float-significand-byte -1)))
98 (defconstant most-positive-short-float most-positive-single-float)
99 (defconstant most-negative-single-float
100 (single-from-bits 1 sb!vm:single-float-normal-exponent-max
101 (ldb sb!vm:single-float-significand-byte -1)))
102 (defconstant most-negative-short-float most-negative-single-float)
103 (defconstant most-positive-double-float
104 (double-from-bits 0 sb!vm:double-float-normal-exponent-max
105 (ldb (byte sb!vm:double-float-digits 0) -1)))
107 (defconstant most-positive-long-float most-positive-double-float)
108 #!+(and long-float x86)
109 (defconstant most-positive-long-float
110 (long-from-bits 0 sb!vm:long-float-normal-exponent-max
111 (ldb (byte sb!vm:long-float-digits 0) -1)))
112 (defconstant most-negative-double-float
113 (double-from-bits 1 sb!vm:double-float-normal-exponent-max
114 (ldb (byte sb!vm:double-float-digits 0) -1)))
116 (defconstant most-negative-long-float most-negative-double-float)
117 #!+(and long-float x86)
118 (defconstant most-negative-long-float
119 (long-from-bits 1 sb!vm:long-float-normal-exponent-max
120 (ldb (byte sb!vm:long-float-digits 0) -1)))
122 ;;; We don't want to do these DEFCONSTANTs at cross-compilation time,
123 ;;; because the cross-compilation host might not support floating
124 ;;; point infinities. Putting them inside a LET removes
125 ;;; toplevel-formness, so that any EVAL-WHEN trickiness in the
126 ;;; DEFCONSTANT forms is suppressed.
128 (defconstant single-float-positive-infinity
129 (single-from-bits 0 (1+ sb!vm:single-float-normal-exponent-max) 0))
130 (defconstant short-float-positive-infinity single-float-positive-infinity)
131 (defconstant single-float-negative-infinity
132 (single-from-bits 1 (1+ sb!vm:single-float-normal-exponent-max) 0))
133 (defconstant short-float-negative-infinity single-float-negative-infinity)
134 (defconstant double-float-positive-infinity
135 (double-from-bits 0 (1+ sb!vm:double-float-normal-exponent-max) 0))
137 (defconstant long-float-positive-infinity double-float-positive-infinity)
138 #!+(and long-float x86)
139 (defconstant long-float-positive-infinity
140 (long-from-bits 0 (1+ sb!vm:long-float-normal-exponent-max)
141 (ash sb!vm:long-float-hidden-bit 32)))
142 (defconstant double-float-negative-infinity
143 (double-from-bits 1 (1+ sb!vm:double-float-normal-exponent-max) 0))
145 (defconstant long-float-negative-infinity double-float-negative-infinity)
146 #!+(and long-float x86)
147 (defconstant long-float-negative-infinity
148 (long-from-bits 1 (1+ sb!vm:long-float-normal-exponent-max)
149 (ash sb!vm:long-float-hidden-bit 32)))
150 ) ; LET-to-suppress-possible-EVAL-WHENs
152 (defconstant single-float-epsilon
153 (single-from-bits 0 (- sb!vm:single-float-bias
154 (1- sb!vm:single-float-digits)) 1))
155 (defconstant short-float-epsilon single-float-epsilon)
156 (defconstant single-float-negative-epsilon
157 (single-from-bits 0 (- sb!vm:single-float-bias sb!vm:single-float-digits) 1))
158 (defconstant short-float-negative-epsilon single-float-negative-epsilon)
159 (defconstant double-float-epsilon
160 (double-from-bits 0 (- sb!vm:double-float-bias
161 (1- sb!vm:double-float-digits)) 1))
163 (defconstant long-float-epsilon double-float-epsilon)
164 #!+(and long-float x86)
165 (defconstant long-float-epsilon
166 (long-from-bits 0 (- sb!vm:long-float-bias (1- sb!vm:long-float-digits))
167 (+ 1 (ash sb!vm:long-float-hidden-bit 32))))
168 (defconstant double-float-negative-epsilon
169 (double-from-bits 0 (- sb!vm:double-float-bias sb!vm:double-float-digits) 1))
171 (defconstant long-float-negative-epsilon double-float-negative-epsilon)
172 #!+(and long-float x86)
173 (defconstant long-float-negative-epsilon
174 (long-from-bits 0 (- sb!vm:long-float-bias sb!vm:long-float-digits)
175 (+ 1 (ash sb!vm:long-float-hidden-bit 32))))
177 ;;;; float predicates and environment query
180 (declaim (maybe-inline float-denormalized-p float-infinity-p float-nan-p
181 float-trapping-nan-p))
183 (defun float-denormalized-p (x)
185 "Return true if the float X is denormalized."
186 (number-dispatch ((x float))
188 (and (zerop (ldb sb!vm:single-float-exponent-byte (single-float-bits x)))
191 (and (zerop (ldb sb!vm:double-float-exponent-byte
192 (double-float-high-bits x)))
194 #!+(and long-float x86)
196 (and (zerop (ldb sb!vm:long-float-exponent-byte (long-float-exp-bits x)))
199 (macrolet ((def (name doc single double #!+(and long-float x86) long)
202 (number-dispatch ((x float))
204 (let ((bits (single-float-bits x)))
205 (and (> (ldb sb!vm:single-float-exponent-byte bits)
206 sb!vm:single-float-normal-exponent-max)
209 (let ((hi (double-float-high-bits x))
210 (lo (double-float-low-bits x)))
211 (declare (ignorable lo))
212 (and (> (ldb sb!vm:double-float-exponent-byte hi)
213 sb!vm:double-float-normal-exponent-max)
215 #!+(and long-float x86)
217 (let ((exp (long-float-exp-bits x))
218 (hi (long-float-high-bits x))
219 (lo (long-float-low-bits x)))
220 (declare (ignorable lo))
221 (and (> (ldb sb!vm:long-float-exponent-byte exp)
222 sb!vm:long-float-normal-exponent-max)
225 (def float-infinity-p
226 "Return true if the float X is an infinity (+ or -)."
227 (zerop (ldb sb!vm:single-float-significand-byte bits))
228 (and (zerop (ldb sb!vm:double-float-significand-byte hi))
230 #!+(and long-float x86)
231 (and (zerop (ldb sb!vm:long-float-significand-byte hi))
235 "Return true if the float X is a NaN (Not a Number)."
236 (not (zerop (ldb sb!vm:single-float-significand-byte bits)))
237 (or (not (zerop (ldb sb!vm:double-float-significand-byte hi)))
239 #!+(and long-float x86)
240 (or (not (zerop (ldb sb!vm:long-float-significand-byte hi)))
243 (def float-trapping-nan-p
244 "Return true if the float X is a trapping NaN (Not a Number)."
245 (zerop (logand (ldb sb!vm:single-float-significand-byte bits)
246 sb!vm:single-float-trapping-nan-bit))
247 (zerop (logand (ldb sb!vm:double-float-significand-byte hi)
248 sb!vm:double-float-trapping-nan-bit))
249 #!+(and long-float x86)
250 (zerop (logand (ldb sb!vm:long-float-significand-byte hi)
251 sb!vm:long-float-trapping-nan-bit))))
253 ;;; If denormalized, use a subfunction from INTEGER-DECODE-FLOAT to find the
254 ;;; actual exponent (and hence how denormalized it is), otherwise we just
255 ;;; return the number of digits or 0.
256 #!-sb-fluid (declaim (maybe-inline float-precision))
257 (defun float-precision (f)
259 "Return a non-negative number of significant digits in its float argument.
260 Will be less than FLOAT-DIGITS if denormalized or zero."
261 (macrolet ((frob (digits bias decode)
263 ((float-denormalized-p f)
264 (multiple-value-bind (ignore exp) (,decode f)
265 (declare (ignore ignore))
267 (+ ,digits (1- ,digits) ,bias exp))))
270 (number-dispatch ((f float))
272 (frob sb!vm:single-float-digits sb!vm:single-float-bias
273 integer-decode-single-denorm))
275 (frob sb!vm:double-float-digits sb!vm:double-float-bias
276 integer-decode-double-denorm))
279 (frob sb!vm:long-float-digits sb!vm:long-float-bias
280 integer-decode-long-denorm)))))
282 (defun float-sign (float1 &optional (float2 (float 1 float1)))
284 "Return a floating-point number that has the same sign as
285 float1 and, if float2 is given, has the same absolute value
287 (declare (float float1 float2))
288 (* (if (etypecase float1
289 (single-float (minusp (single-float-bits float1)))
290 (double-float (minusp (double-float-high-bits float1)))
292 (long-float (minusp (long-float-exp-bits float1))))
297 (defun float-format-digits (format)
299 ((short-float single-float) sb!vm:single-float-digits)
300 ((double-float #!-long-float long-float) sb!vm:double-float-digits)
302 (long-float sb!vm:long-float-digits)))
304 #!-sb-fluid (declaim (inline float-digits float-radix))
306 (defun float-digits (f)
307 (number-dispatch ((f float))
308 ((single-float) sb!vm:single-float-digits)
309 ((double-float) sb!vm:double-float-digits)
311 ((long-float) sb!vm:long-float-digits)))
313 (setf (fdefinition 'float-radix)
314 ;; FIXME: Python flushes unused variable X in CLAMBDA, then
315 ;; flushes unused reference to X in XEP together with type
316 ;; check. When this is fixed, rewrite this definition in an
317 ;; ordinary form. -- APD, 2002-10-21
320 "Return (as an integer) the radix b of its floating-point argument."
322 (error 'type-error :datum x :expected-type 'float))
325 ;;;; INTEGER-DECODE-FLOAT and DECODE-FLOAT
328 (declaim (maybe-inline integer-decode-single-float
329 integer-decode-double-float))
331 ;;; Handle the denormalized case of INTEGER-DECODE-FLOAT for SINGLE-FLOAT.
332 (defun integer-decode-single-denorm (x)
333 (declare (type single-float x))
334 (let* ((bits (single-float-bits (abs x)))
335 (sig (ash (ldb sb!vm:single-float-significand-byte bits) 1))
337 (declare (type (unsigned-byte 24) sig)
338 (type (integer 0 23) extra-bias))
340 (unless (zerop (logand sig sb!vm:single-float-hidden-bit))
342 (setq sig (ash sig 1))
345 (- (- sb!vm:single-float-bias)
346 sb!vm:single-float-digits
348 (if (minusp (float-sign x)) -1 1))))
350 ;;; Handle the single-float case of INTEGER-DECODE-FLOAT. If an infinity or
351 ;;; NaN, error. If a denorm, call i-d-s-DENORM to handle it.
352 (defun integer-decode-single-float (x)
353 (declare (single-float x))
354 (let* ((bits (single-float-bits (abs x)))
355 (exp (ldb sb!vm:single-float-exponent-byte bits))
356 (sig (ldb sb!vm:single-float-significand-byte bits))
357 (sign (if (minusp (float-sign x)) -1 1))
358 (biased (- exp sb!vm:single-float-bias sb!vm:single-float-digits)))
359 (declare (fixnum biased))
360 (unless (<= exp sb!vm:single-float-normal-exponent-max)
361 (error "can't decode NaN or infinity: ~S" x))
362 (cond ((and (zerop exp) (zerop sig))
363 (values 0 biased sign))
364 ((< exp sb!vm:single-float-normal-exponent-min)
365 (integer-decode-single-denorm x))
367 (values (logior sig sb!vm:single-float-hidden-bit) biased sign)))))
369 ;;; like INTEGER-DECODE-SINGLE-DENORM, only doubly so
370 (defun integer-decode-double-denorm (x)
371 (declare (type double-float x))
372 (let* ((high-bits (double-float-high-bits (abs x)))
373 (sig-high (ldb sb!vm:double-float-significand-byte high-bits))
374 (low-bits (double-float-low-bits x))
375 (sign (if (minusp (float-sign x)) -1 1))
376 (biased (- (- sb!vm:double-float-bias) sb!vm:double-float-digits)))
379 (extra-bias (- sb!vm:double-float-digits 33))
381 (declare (type (unsigned-byte 32) sig) (fixnum extra-bias))
383 (unless (zerop (logand sig bit)) (return))
384 (setq sig (ash sig 1))
386 (values (ash sig (- sb!vm:double-float-digits 32))
387 (truly-the fixnum (- biased extra-bias))
389 (let ((sig (ash sig-high 1))
391 (declare (type (unsigned-byte 32) sig) (fixnum extra-bias))
393 (unless (zerop (logand sig sb!vm:double-float-hidden-bit))
395 (setq sig (ash sig 1))
397 (values (logior (ash sig 32) (ash low-bits (1- extra-bias)))
398 (truly-the fixnum (- biased extra-bias))
401 ;;; like INTEGER-DECODE-SINGLE-FLOAT, only doubly so
402 (defun integer-decode-double-float (x)
403 (declare (double-float x))
405 (hi (double-float-high-bits abs))
406 (lo (double-float-low-bits abs))
407 (exp (ldb sb!vm:double-float-exponent-byte hi))
408 (sig (ldb sb!vm:double-float-significand-byte hi))
409 (sign (if (minusp (float-sign x)) -1 1))
410 (biased (- exp sb!vm:double-float-bias sb!vm:double-float-digits)))
411 (declare (fixnum biased))
412 (unless (<= exp sb!vm:double-float-normal-exponent-max)
413 (error "Can't decode NaN or infinity: ~S." x))
414 (cond ((and (zerop exp) (zerop sig) (zerop lo))
415 (values 0 biased sign))
416 ((< exp sb!vm:double-float-normal-exponent-min)
417 (integer-decode-double-denorm x))
420 (logior (ash (logior (ldb sb!vm:double-float-significand-byte hi)
421 sb!vm:double-float-hidden-bit)
426 #!+(and long-float x86)
427 (defun integer-decode-long-denorm (x)
428 (declare (type long-float x))
429 (let* ((high-bits (long-float-high-bits (abs x)))
430 (sig-high (ldb sb!vm:long-float-significand-byte high-bits))
431 (low-bits (long-float-low-bits x))
432 (sign (if (minusp (float-sign x)) -1 1))
433 (biased (- (- sb!vm:long-float-bias) sb!vm:long-float-digits)))
436 (extra-bias (- sb!vm:long-float-digits 33))
438 (declare (type (unsigned-byte 32) sig) (fixnum extra-bias))
440 (unless (zerop (logand sig bit)) (return))
441 (setq sig (ash sig 1))
443 (values (ash sig (- sb!vm:long-float-digits 32))
444 (truly-the fixnum (- biased extra-bias))
446 (let ((sig (ash sig-high 1))
448 (declare (type (unsigned-byte 32) sig) (fixnum extra-bias))
450 (unless (zerop (logand sig sb!vm:long-float-hidden-bit))
452 (setq sig (ash sig 1))
454 (values (logior (ash sig 32) (ash low-bits (1- extra-bias)))
455 (truly-the fixnum (- biased extra-bias))
458 #!+(and long-float x86)
459 (defun integer-decode-long-float (x)
460 (declare (long-float x))
461 (let* ((hi (long-float-high-bits x))
462 (lo (long-float-low-bits x))
463 (exp-bits (long-float-exp-bits x))
464 (exp (ldb sb!vm:long-float-exponent-byte exp-bits))
465 (sign (if (minusp exp-bits) -1 1))
466 (biased (- exp sb!vm:long-float-bias sb!vm:long-float-digits)))
467 (declare (fixnum biased))
468 (unless (<= exp sb!vm:long-float-normal-exponent-max)
469 (error "can't decode NaN or infinity: ~S" x))
470 (cond ((and (zerop exp) (zerop hi) (zerop lo))
471 (values 0 biased sign))
472 ((< exp sb!vm:long-float-normal-exponent-min)
473 (integer-decode-long-denorm x))
475 (values (logior (ash hi 32) lo) biased sign)))))
477 ;;; Dispatch to the correct type-specific i-d-f function.
478 (defun integer-decode-float (x)
480 "Return three values:
481 1) an integer representation of the significand.
482 2) the exponent for the power of 2 that the significand must be multiplied
483 by to get the actual value. This differs from the DECODE-FLOAT exponent
484 by FLOAT-DIGITS, since the significand has been scaled to have all its
485 digits before the radix point.
486 3) -1 or 1 (i.e. the sign of the argument.)"
487 (number-dispatch ((x float))
489 (integer-decode-single-float x))
491 (integer-decode-double-float x))
494 (integer-decode-long-float x))))
496 #!-sb-fluid (declaim (maybe-inline decode-single-float decode-double-float))
498 ;;; Handle the denormalized case of DECODE-SINGLE-FLOAT. We call
499 ;;; INTEGER-DECODE-SINGLE-DENORM and then make the result into a float.
500 (defun decode-single-denorm (x)
501 (declare (type single-float x))
502 (multiple-value-bind (sig exp sign) (integer-decode-single-denorm x)
503 (values (make-single-float
504 (dpb sig sb!vm:single-float-significand-byte
505 (dpb sb!vm:single-float-bias
506 sb!vm:single-float-exponent-byte
508 (truly-the fixnum (+ exp sb!vm:single-float-digits))
511 ;;; Handle the single-float case of DECODE-FLOAT. If an infinity or NaN,
512 ;;; error. If a denorm, call d-s-DENORM to handle it.
513 (defun decode-single-float (x)
514 (declare (single-float x))
515 (let* ((bits (single-float-bits (abs x)))
516 (exp (ldb sb!vm:single-float-exponent-byte bits))
517 (sign (float-sign x))
518 (biased (truly-the single-float-exponent
519 (- exp sb!vm:single-float-bias))))
520 (unless (<= exp sb!vm:single-float-normal-exponent-max)
521 (error "can't decode NaN or infinity: ~S" x))
523 (values 0.0f0 biased sign))
524 ((< exp sb!vm:single-float-normal-exponent-min)
525 (decode-single-denorm x))
527 (values (make-single-float
528 (dpb sb!vm:single-float-bias
529 sb!vm:single-float-exponent-byte
533 ;;; like DECODE-SINGLE-DENORM, only doubly so
534 (defun decode-double-denorm (x)
535 (declare (double-float x))
536 (multiple-value-bind (sig exp sign) (integer-decode-double-denorm x)
537 (values (make-double-float
538 (dpb (logand (ash sig -32) (lognot sb!vm:double-float-hidden-bit))
539 sb!vm:double-float-significand-byte
540 (dpb sb!vm:double-float-bias
541 sb!vm:double-float-exponent-byte 0))
542 (ldb (byte 32 0) sig))
543 (truly-the fixnum (+ exp sb!vm:double-float-digits))
546 ;;; like DECODE-SINGLE-FLOAT, only doubly so
547 (defun decode-double-float (x)
548 (declare (double-float x))
550 (hi (double-float-high-bits abs))
551 (lo (double-float-low-bits abs))
552 (exp (ldb sb!vm:double-float-exponent-byte hi))
553 (sign (float-sign x))
554 (biased (truly-the double-float-exponent
555 (- exp sb!vm:double-float-bias))))
556 (unless (<= exp sb!vm:double-float-normal-exponent-max)
557 (error "can't decode NaN or infinity: ~S" x))
559 (values 0.0d0 biased sign))
560 ((< exp sb!vm:double-float-normal-exponent-min)
561 (decode-double-denorm x))
563 (values (make-double-float
564 (dpb sb!vm:double-float-bias
565 sb!vm:double-float-exponent-byte hi)
569 #!+(and long-float x86)
570 (defun decode-long-denorm (x)
571 (declare (long-float x))
572 (multiple-value-bind (sig exp sign) (integer-decode-long-denorm x)
573 (values (make-long-float sb!vm:long-float-bias (ash sig -32)
574 (ldb (byte 32 0) sig))
575 (truly-the fixnum (+ exp sb!vm:long-float-digits))
578 #!+(and long-float x86)
579 (defun decode-long-float (x)
580 (declare (long-float x))
581 (let* ((hi (long-float-high-bits x))
582 (lo (long-float-low-bits x))
583 (exp-bits (long-float-exp-bits x))
584 (exp (ldb sb!vm:long-float-exponent-byte exp-bits))
585 (sign (if (minusp exp-bits) -1l0 1l0))
586 (biased (truly-the long-float-exponent
587 (- exp sb!vm:long-float-bias))))
588 (unless (<= exp sb!vm:long-float-normal-exponent-max)
589 (error "can't decode NaN or infinity: ~S" x))
591 (values 0.0l0 biased sign))
592 ((< exp sb!vm:long-float-normal-exponent-min)
593 (decode-long-denorm x))
595 (values (make-long-float
596 (dpb sb!vm:long-float-bias sb!vm:long-float-exponent-byte
602 ;;; Dispatch to the appropriate type-specific function.
603 (defun decode-float (f)
605 "Return three values:
606 1) a floating-point number representing the significand. This is always
607 between 0.5 (inclusive) and 1.0 (exclusive).
608 2) an integer representing the exponent.
609 3) -1.0 or 1.0 (i.e. the sign of the argument.)"
610 (number-dispatch ((f float))
612 (decode-single-float f))
614 (decode-double-float f))
617 (decode-long-float f))))
621 #!-sb-fluid (declaim (maybe-inline scale-single-float scale-double-float))
623 ;;; Handle float scaling where the X is denormalized or the result is
624 ;;; denormalized or underflows to 0.
625 (defun scale-float-maybe-underflow (x exp)
626 (multiple-value-bind (sig old-exp) (integer-decode-float x)
627 (let* ((digits (float-digits x))
628 (new-exp (+ exp old-exp digits
630 (single-float sb!vm:single-float-bias)
631 (double-float sb!vm:double-float-bias))))
632 (sign (if (minusp (float-sign x)) 1 0)))
636 (single-float sb!vm:single-float-normal-exponent-min)
637 (double-float sb!vm:double-float-normal-exponent-min)))
638 (when (sb!vm:current-float-trap :inexact)
639 (error 'floating-point-inexact :operation 'scale-float
640 :operands (list x exp)))
641 (when (sb!vm:current-float-trap :underflow)
642 (error 'floating-point-underflow :operation 'scale-float
643 :operands (list x exp)))
644 (let ((shift (1- new-exp)))
645 (if (< shift (- (1- digits)))
648 (single-float (single-from-bits sign 0 (ash sig shift)))
649 (double-float (double-from-bits sign 0 (ash sig shift)))))))
652 (single-float (single-from-bits sign new-exp sig))
653 (double-float (double-from-bits sign new-exp sig))))))))
655 ;;; Called when scaling a float overflows, or the original float was a
656 ;;; NaN or infinity. If overflow errors are trapped, then error,
657 ;;; otherwise return the appropriate infinity. If a NaN, signal or not
659 (defun scale-float-maybe-overflow (x exp)
661 ((float-infinity-p x)
662 ;; Infinity is infinity, no matter how small...
665 (when (and (float-trapping-nan-p x)
666 (sb!vm:current-float-trap :invalid))
667 (error 'floating-point-invalid-operation :operation 'scale-float
668 :operands (list x exp)))
671 (when (sb!vm:current-float-trap :overflow)
672 (error 'floating-point-overflow :operation 'scale-float
673 :operands (list x exp)))
674 (when (sb!vm:current-float-trap :inexact)
675 (error 'floating-point-inexact :operation 'scale-float
676 :operands (list x exp)))
679 (single-float single-float-positive-infinity)
680 (double-float double-float-positive-infinity))))))
682 ;;; Scale a single or double float, calling the correct over/underflow
684 (defun scale-single-float (x exp)
685 (declare (single-float x) (fixnum exp))
686 (let* ((bits (single-float-bits x))
687 (old-exp (ldb sb!vm:single-float-exponent-byte bits))
688 (new-exp (+ old-exp exp)))
691 ((or (< old-exp sb!vm:single-float-normal-exponent-min)
692 (< new-exp sb!vm:single-float-normal-exponent-min))
693 (scale-float-maybe-underflow x exp))
694 ((or (> old-exp sb!vm:single-float-normal-exponent-max)
695 (> new-exp sb!vm:single-float-normal-exponent-max))
696 (scale-float-maybe-overflow x exp))
698 (make-single-float (dpb new-exp
699 sb!vm:single-float-exponent-byte
701 (defun scale-double-float (x exp)
702 (declare (double-float x) (fixnum exp))
703 (let* ((hi (double-float-high-bits x))
704 (lo (double-float-low-bits x))
705 (old-exp (ldb sb!vm:double-float-exponent-byte hi))
706 (new-exp (+ old-exp exp)))
709 ((or (< old-exp sb!vm:double-float-normal-exponent-min)
710 (< new-exp sb!vm:double-float-normal-exponent-min))
711 (scale-float-maybe-underflow x exp))
712 ((or (> old-exp sb!vm:double-float-normal-exponent-max)
713 (> new-exp sb!vm:double-float-normal-exponent-max))
714 (scale-float-maybe-overflow x exp))
716 (make-double-float (dpb new-exp sb!vm:double-float-exponent-byte hi)
719 #!+(and x86 long-float)
720 (defun scale-long-float (x exp)
721 (declare (long-float x) (fixnum exp))
724 ;;; Dispatch to the correct type-specific scale-float function.
725 (defun scale-float (f ex)
727 "Return the value (* f (expt (float 2 f) ex)), but with no unnecessary loss
728 of precision or overflow."
729 (number-dispatch ((f float))
731 (scale-single-float f ex))
733 (scale-double-float f ex))
736 (scale-long-float f ex))))
738 ;;;; converting to/from floats
740 (defun float (number &optional (other () otherp))
742 "Converts any REAL to a float. If OTHER is not provided, it returns a
743 SINGLE-FLOAT if NUMBER is not already a FLOAT. If OTHER is provided, the
744 result is the same float format as OTHER."
746 (number-dispatch ((number real) (other float))
747 (((foreach rational single-float double-float #!+long-float long-float)
748 (foreach single-float double-float #!+long-float long-float))
749 (coerce number '(dispatch-type other))))
752 (coerce number 'single-float))))
754 (macrolet ((frob (name type)
756 (number-dispatch ((x real))
757 (((foreach single-float double-float #!+long-float long-float
761 (bignum-to-float x ',type))
763 (float-ratio x ',type))))))
764 (frob %single-float single-float)
765 (frob %double-float double-float)
767 (frob %long-float long-float))
769 ;;; Convert a ratio to a float. We avoid any rounding error by doing an
770 ;;; integer division. Accuracy is important to preserve read/print
771 ;;; consistency, since this is ultimately how the reader reads a float. We
772 ;;; scale the numerator by a power of two until the division results in the
773 ;;; desired number of fraction bits, then do round-to-nearest.
774 (defun float-ratio (x format)
775 (let* ((signed-num (numerator x))
776 (plusp (plusp signed-num))
777 (num (if plusp signed-num (- signed-num)))
778 (den (denominator x))
779 (digits (float-format-digits format))
781 (declare (fixnum digits scale))
782 ;; Strip any trailing zeros from the denominator and move it into the scale
783 ;; factor (to minimize the size of the operands.)
784 (let ((den-twos (1- (integer-length (logxor den (1- den))))))
785 (declare (fixnum den-twos))
786 (decf scale den-twos)
787 (setq den (ash den (- den-twos))))
788 ;; Guess how much we need to scale by from the magnitudes of the numerator
789 ;; and denominator. We want one extra bit for a guard bit.
790 (let* ((num-len (integer-length num))
791 (den-len (integer-length den))
792 (delta (- den-len num-len))
793 (shift (1+ (the fixnum (+ delta digits))))
794 (shifted-num (ash num shift)))
795 (declare (fixnum delta shift))
797 (labels ((float-and-scale (bits)
798 (let* ((bits (ash bits -1))
799 (len (integer-length bits)))
800 (cond ((> len digits)
801 (aver (= len (the fixnum (1+ digits))))
802 (scale-float (floatit (ash bits -1)) (1+ scale)))
804 (scale-float (floatit bits) scale)))))
806 (let ((sign (if plusp 0 1)))
809 (single-from-bits sign sb!vm:single-float-bias bits))
811 (double-from-bits sign sb!vm:double-float-bias bits))
814 (long-from-bits sign sb!vm:long-float-bias bits))))))
816 (multiple-value-bind (fraction-and-guard rem)
817 (truncate shifted-num den)
818 (let ((extra (- (integer-length fraction-and-guard) digits)))
819 (declare (fixnum extra))
822 ((oddp fraction-and-guard)
826 (if (zerop (logand fraction-and-guard 2))
828 (1+ fraction-and-guard)))
829 (float-and-scale (1+ fraction-and-guard)))))
831 (return (float-and-scale fraction-and-guard)))))
832 (setq shifted-num (ash shifted-num -1))
836 These might be useful if we ever have a machine without float/integer
837 conversion hardware. For now, we'll use special ops that
838 uninterruptibly frob the rounding modes & do ieee round-to-integer.
840 ;;; The compiler compiles a call to this when we are doing %UNARY-TRUNCATE
841 ;;; and the result is known to be a fixnum. We can avoid some generic
842 ;;; arithmetic in this case.
843 (defun %unary-truncate-single-float/fixnum (x)
844 (declare (single-float x) (values fixnum))
845 (locally (declare (optimize (speed 3) (safety 0)))
846 (let* ((bits (single-float-bits x))
847 (exp (ldb sb!vm:single-float-exponent-byte bits))
848 (frac (logior (ldb sb!vm:single-float-significand-byte bits)
849 sb!vm:single-float-hidden-bit))
850 (shift (- exp sb!vm:single-float-digits sb!vm:single-float-bias)))
851 (when (> exp sb!vm:single-float-normal-exponent-max)
852 (error 'floating-point-invalid-operation :operator 'truncate
854 (if (<= shift (- sb!vm:single-float-digits))
856 (let ((res (ash frac shift)))
857 (declare (type (unsigned-byte 31) res))
862 ;;; Double-float version of this operation (see above single op).
863 (defun %unary-truncate-double-float/fixnum (x)
864 (declare (double-float x) (values fixnum))
865 (locally (declare (optimize (speed 3) (safety 0)))
866 (let* ((hi-bits (double-float-high-bits x))
867 (exp (ldb sb!vm:double-float-exponent-byte hi-bits))
868 (frac (logior (ldb sb!vm:double-float-significand-byte hi-bits)
869 sb!vm:double-float-hidden-bit))
870 (shift (- exp (- sb!vm:double-float-digits sb!vm:n-word-bits)
871 sb!vm:double-float-bias)))
872 (when (> exp sb!vm:double-float-normal-exponent-max)
873 (error 'floating-point-invalid-operation :operator 'truncate
875 (if (<= shift (- sb!vm:n-word-bits sb!vm:double-float-digits))
877 (let* ((res-hi (ash frac shift))
878 (res (if (plusp shift)
881 (ash (double-float-low-bits x)
882 (- shift sb!vm:n-word-bits))))
884 (declare (type (unsigned-byte 31) res-hi res))
890 ;;; This function is called when we are doing a truncate without any funky
891 ;;; divisor, i.e. converting a float or ratio to an integer. Note that we do
892 ;;; *not* return the second value of truncate, so it must be computed by the
893 ;;; caller if needed.
895 ;;; In the float case, we pick off small arguments so that compiler can use
896 ;;; special-case operations. We use an exclusive test, since (due to round-off
897 ;;; error), (float most-positive-fixnum) may be greater than
898 ;;; most-positive-fixnum.
899 (defun %unary-truncate (number)
900 (number-dispatch ((number real))
902 ((ratio) (values (truncate (numerator number) (denominator number))))
903 (((foreach single-float double-float #!+long-float long-float))
904 (if (< (float most-negative-fixnum number)
906 (float most-positive-fixnum number))
907 (truly-the fixnum (%unary-truncate number))
908 (multiple-value-bind (bits exp) (integer-decode-float number)
909 (let ((res (ash bits exp)))
914 ;;; Similar to %UNARY-TRUNCATE, but rounds to the nearest integer. If we
915 ;;; can't use the round primitive, then we do our own round-to-nearest on the
916 ;;; result of i-d-f. [Note that this rounding will really only happen with
917 ;;; double floats, since the whole single-float fraction will fit in a fixnum,
918 ;;; so all single-floats larger than most-positive-fixnum can be precisely
919 ;;; represented by an integer.]
920 (defun %unary-round (number)
921 (number-dispatch ((number real))
923 ((ratio) (values (round (numerator number) (denominator number))))
924 (((foreach single-float double-float #!+long-float long-float))
925 (if (< (float most-negative-fixnum number)
927 (float most-positive-fixnum number))
928 (truly-the fixnum (%unary-round number))
929 (multiple-value-bind (bits exp) (integer-decode-float number)
930 (let* ((shifted (ash bits exp))
931 (rounded (if (and (minusp exp)
934 (lognot (ash -1 (- exp))))
944 "RATIONAL produces a rational number for any real numeric argument. This is
945 more efficient than RATIONALIZE, but it assumes that floating-point is
946 completely accurate, giving a result that isn't as pretty."
947 (number-dispatch ((x real))
948 (((foreach single-float double-float #!+long-float long-float))
949 (multiple-value-bind (bits exp) (integer-decode-float x)
952 (let* ((int (if (minusp x) (- bits) bits))
953 (digits (float-digits x))
956 (integer-/-integer int (ash 1 (+ digits (- ex))))
957 (integer-/-integer (ash int ex) (ash 1 digits)))))))
960 (defun rationalize (x)
962 "Converts any REAL to a RATIONAL. Floats are converted to a simple rational
963 representation exploiting the assumption that floats are only accurate to
964 their precision. RATIONALIZE (and also RATIONAL) preserve the invariant:
965 (= x (float (rationalize x) x))"
966 (number-dispatch ((x real))
967 (((foreach single-float double-float #!+long-float long-float))
968 ;; Thanks to Kim Fateman, who stole this function rationalize-float from
969 ;; macsyma's rational. Macsyma'a rationalize was written by the legendary
970 ;; Gosper (rwg). Guy Steele said about Gosper, "He has been called the
971 ;; only living 17th century mathematician and is also the best pdp-10
972 ;; hacker I know." So, if you can understand or debug this code you win
974 (cond ((minusp x) (- (rationalize (- x))))
977 (let ((eps (etypecase x
978 (single-float single-float-epsilon)
979 (double-float double-float-epsilon)
981 (long-float long-float-epsilon)))
984 (do ((xx x (setq y (/ (float 1.0 x) (- xx (float a x)))))
985 (num (setq a (truncate x))
986 (+ (* (setq a (truncate y)) num) onum))
987 (den 1 (+ (* a den) oden))
990 ((and (not (zerop den))
991 (not (> (abs (/ (- x (/ (float num x)
995 (integer-/-integer num den))
996 (declare ((dispatch-type x) xx)))))))