# Diff of /src/code/irrat.lisp

revision 1.41 by rtoy, Tue Oct 19 15:07:30 2004 UTC revision 1.41.2.1 by rtoy, Mon Dec 19 01:09:51 2005 UTC
# Line 373  Line 373
373                        x)))                        x)))
374            (+ n (log (scale-float (float f 1d0) (- exp))            (+ n (log (scale-float (float f 1d0) (- exp))
375                      2d0))))))                      2d0))))))
376
377  (defun log (number &optional (base nil base-p))  (defun log (number &optional (base nil base-p))
378    "Return the logarithm of NUMBER in the base BASE, which defaults to e."    "Return the logarithm of NUMBER in the base BASE, which defaults to e."
379    (if base-p    (if base-p
# Line 386  Line 386
386               ;; and the base are positive integers.  Use the rule that               ;; and the base are positive integers.  Use the rule that
387               ;; log_b(x) = log_2(x)/log_2(b)               ;; log_b(x) = log_2(x)/log_2(b)
388               (coerce (/ (log2 number) (log2 base)) 'single-float))               (coerce (/ (log2 number) (log2 base)) 'single-float))
389                ((and (realp number) (realp base))
390                 ;; CLHS 12.1.4.1 says
391                 ;;
392                 ;;   When rationals and floats are combined by a
393                 ;;   numerical function, the rational is first converted
394                 ;;   to a float of the same format.
395                 ;;
396                 ;; So assume this applies to floats as well convert all
397                 ;; numbers to the largest float format before computing
398                 ;; the log.
399                 ;;
400                 ;; This makes (log 17 10.0) = (log 17.0 10) and so on.
401                 (number-dispatch ((number real) (base real))
402                   ((double-float
403                     (foreach double-float single-float fixnum bignum ratio))
404                    (/ (log number) (log (coerce base 'double-float))))
405                   (((foreach single-float fixnum bignum ratio)
406                     double-float)
407                    (/ (log (coerce number 'double-float)) (log base)))
408                   (((foreach single-float fixnum bignum ratio)
409                     (foreach single-float fixnum bignum ratio))
410                    ;; Converting everything to double-float helps the
411                    ;; cases like (log 17 10) = (/ (log 17) (log 10)).
412                    ;; This is usually handled above, but if we compute (/
413                    ;; (log 17) (log 10)), we get a slightly different
414                    ;; answer due to roundoff.  This makes it a bit more
415                    ;; consistent.
416                    ;;
417                    ;; FIXME: This probably needs more work.
418                    (let ((result (/ (log (float number 1d0))
419                                     (log (float base 1d0)))))
420                      (if (realp result)
421                          (coerce result 'single-float)
422                          (coerce result '(complex single-float)))))))
423              (t              (t
424                 ;; FIXME:  This probably needs some work as well.
425               (/ (log number) (log base))))               (/ (log number) (log base))))
426        (number-dispatch ((number number))        (number-dispatch ((number number))
427          (((foreach fixnum bignum))          (((foreach fixnum bignum))
# Line 531  Line 566
566           (complex-asin number)           (complex-asin number)
567           (coerce (%asin (coerce number 'double-float)) 'single-float)))           (coerce (%asin (coerce number 'double-float)) 'single-float)))
568      (((foreach single-float double-float))      (((foreach single-float double-float))
569       (if (or (> number (coerce 1 '(dispatch-type number)))       (if (or (float-nan-p number)
570               (< number (coerce -1 '(dispatch-type number))))               (and (<= number (coerce 1 '(dispatch-type number)))
571           (complex-asin number)                    (>= number (coerce -1 '(dispatch-type number)))))
572           (coerce (%asin (coerce number 'double-float))           (coerce (%asin (coerce number 'double-float))
573                   '(dispatch-type number))))                   '(dispatch-type number))
574                     (complex-asin number)))
575      ((complex)      ((complex)
576       (complex-asin number))))       (complex-asin number))))
577
# Line 547  Line 583
583           (complex-acos number)           (complex-acos number)
584           (coerce (%acos (coerce number 'double-float)) 'single-float)))           (coerce (%acos (coerce number 'double-float)) 'single-float)))
585      (((foreach single-float double-float))      (((foreach single-float double-float))
586       (if (or (> number (coerce 1 '(dispatch-type number)))       (if (or (float-nan-p number)
587               (< number (coerce -1 '(dispatch-type number))))               (and (<= number (coerce 1 '(dispatch-type number)))
588           (complex-acos number)                    (>= number (coerce -1 '(dispatch-type number)))))
589           (coerce (%acos (coerce number 'double-float))           (coerce (%acos (coerce number 'double-float))
590                   '(dispatch-type number))))                   '(dispatch-type number))
591                     (complex-acos number)))
592      ((complex)      ((complex)
593       (complex-acos number))))       (complex-acos number))))
594

Legend:
 Removed from v.1.41 changed lines Added in v.1.41.2.1