# Diff of /src/code/irrat.lisp

revision 1.55 by rtoy, Thu Jan 31 19:12:40 2008 UTC revision 1.56 by rtoy, Wed Oct 8 17:03:07 2008 UTC
# Line 455  Line 455
455                                     (coerce (* pow (%cos y*pi)) rtype)                                     (coerce (* pow (%cos y*pi)) rtype)
456                                     (coerce (* pow (%sin y*pi)) rtype)))))))))))))                                     (coerce (* pow (%sin y*pi)) rtype)))))))))))))
457        (declare (inline real-expt))        (declare (inline real-expt))
458          ;; This is really messy and should be cleaned up.  The easiest
459          ;; way to see if we're doing what we should is the macroexpand
460          ;; the number-dispatch and check each branch.
461          ;;
462          ;; We try to apply the rule of float precision contagion (CLHS
463          ;; 12.1.4.4): the result has the same precision has the most
464          ;; precise argument.
465        (number-dispatch ((base number) (power number))        (number-dispatch ((base number) (power number))
466          (((foreach fixnum (or bignum ratio) (complex rational)) integer)          (((foreach fixnum (or bignum ratio) (complex rational))
467              integer)
468           (intexp base power))           (intexp base power))
469          (((foreach single-float double-float) rational)          (((foreach single-float double-float)
470              rational)
471           (real-expt base power '(dispatch-type base)))           (real-expt base power '(dispatch-type base)))
472          (((foreach fixnum (or bignum ratio) single-float)          (((foreach fixnum (or bignum ratio) single-float)
473            (foreach ratio single-float))            (foreach ratio single-float))
# Line 469  Line 478
478          ((double-float single-float)          ((double-float single-float)
479           (real-expt base power 'double-float))           (real-expt base power 'double-float))
480          #+double-double          #+double-double
481          (((foreach fixnum (or bignum ratio) single-float double-float double-double-float)          (((foreach fixnum (or bignum ratio) single-float double-float
482                       double-double-float)
483            double-double-float)            double-double-float)
484           (dd-%pow (coerce base 'double-double-float) power))           (dd-%pow (coerce base 'double-double-float) power))
485          #+double-double          #+double-double
486          ((double-double-float          ((double-double-float
487            (foreach fixnum (or bignum ratio) single-float double-float))            (foreach fixnum (or bignum ratio) single-float double-float))
488           (dd-%pow base (coerce power 'double-double-float)))           (dd-%pow base (coerce power 'double-double-float)))
489          (((foreach (complex rational) (complex float)) rational)          (((foreach (complex rational) (complex single-float) (complex double-float)
490                       #+double-double (complex double-double-float))
491              rational)
492           (* (expt (abs base) power)           (* (expt (abs base) power)
493              (cis (* power (phase base)))))              (cis (* power (phase base)))))
494          (((foreach fixnum (or bignum ratio) single-float double-float          #+double-double
495                     #+double-double double-double-float)          ((double-double-float
496            complex)            complex)
497           (if (and (zerop base) (plusp (realpart power)))           (if (and (zerop base) (plusp (realpart power)))
498               (* base power)               (* base power)
499                 (exp (* power (* (log2 base 1w0) (log 2w0))))))
500            (((foreach fixnum (or bignum ratio) single-float double-float)
501              (foreach (complex double-float)))
502             ;; Result should have double-float accuracy.  Use log2 in
503             ;; case the base won't fit in a double-float.
504             (if (and (zerop base) (plusp (realpart power)))
505                 (* base power)
506                 (exp (* power (* (log2 base) (log 2d0))))))
507            ((double-float
508              (foreach (complex rational) (complex single-float)))
509             (if (and (zerop base) (plusp (realpart power)))
510                 (* base power)
511                 (exp (* power (log base)))))
512            #+double-double
513            (((foreach fixnum (or bignum ratio) single-float double-float)
514              (foreach (complex double-double-float)))
515             ;; Result should have double-double-float accuracy.  Use log2
516             ;; in case the base won't fit in a double-float.
517             (if (and (zerop base) (plusp (realpart power)))
518                 (* base power)
519                 (exp (* power (* (log2 base 1w0) (log 2w0))))))
520            (((foreach fixnum (or bignum ratio) single-float)
521              (foreach (complex single-float)))
522             (if (and (zerop base) (plusp (realpart power)))
523                 (* base power)
524                 (exp (* power (log base)))))
525            (((foreach (complex rational) (complex single-float))
526              (foreach single-float (complex single-float)))
527             (if (and (zerop base) (plusp (realpart power)))
528                 (* base power)
529                 (exp (* power (log base)))))
530            (((foreach (complex rational) (complex single-float))
531              (foreach double-float (complex double-float)))
532             (if (and (zerop base) (plusp (realpart power)))
533                 (* base power)
534                 (exp (* power (log (coerce base '(complex double-float)))))))
535            #+double-double
536            (((foreach (complex rational) (complex single-float))
537              (foreach double-double-float (complex double-double-float)))
538             (if (and (zerop base) (plusp (realpart power)))
539                 (* base power)
540                 (exp (* power (log (coerce base '(complex double-double-float)))))))
541            (((foreach (complex double-float))
542              (foreach single-float double-float (complex single-float)
543                       (complex double-float)))
544             (if (and (zerop base) (plusp (realpart power)))
545                 (* base power)
546               (exp (* power (log base)))))               (exp (* power (log base)))))
547          (((foreach (complex float) (complex rational))          #+double-double
548            (foreach complex double-float single-float #+double-double double-double-float))          (((foreach (complex double-float))
549              (foreach double-double-float (complex double-double-float)))
550             (if (and (zerop base) (plusp (realpart power)))
551                 (* base power)
552                 (exp (* power (log (coerce base '(complex double-double-float)))))))
553            #+double-double
554            (((foreach (complex double-double-float))
555              (foreach float (complex float)))
556           (if (and (zerop base) (plusp (realpart power)))           (if (and (zerop base) (plusp (realpart power)))
557               (* base power)               (* base power)
558               (exp (* power (log base)))))))))               (exp (* power (log base)))))))))

Legend:
 Removed from v.1.55 changed lines Added in v.1.56