re-align expt' and flexpt' to match C99 pow() spec

Also, improve precision of some complex results to avoid
excessive `+nan.0's.

Closes PR 12935
This commit is contained in:
Matthew Flatt 2012-08-06 11:07:04 -06:00
parent 7494fccc4c
commit 13d7a37eb6
12 changed files with 460 additions and 109 deletions

View File

@ -80,9 +80,31 @@ or @racket[flsqrt].}
flonum?]{ flonum?]{
Like @racket[expt], but constrained to consume and produce Like @racket[expt], but constrained to consume and produce
@tech{flonums}. The result is @racket[+nan.0] when @racket[a] is @tech{flonums}.
negative and @racket[b] is not an integer or when @racket[a] is zero
and @racket[b] is not positive.} Due to the result constraint, the results compared to @racket[expt]
differ in the following cases:
@margin-note*{These special cases correspond to @tt{pow} in C99 @cite["C99"].}
@;
@itemlist[#:style 'compact
@item{@racket[(flexpt -1.0 +inf.0)] --- @racket[1.0]}
@item{@racket[(flexpt a +inf.0)] where @racket[a] is
negative --- @racket[(expt (abs a) +inf.0)]}
@item{@racket[(flexpt a -inf.0)] where @racket[a] is
negative --- @racket[(expt (abs a) -inf.0)]}
@item{@racket[(expt -inf.0 b)] where @racket[b] is a non-integer:
@itemlist[#:style 'compact
@item{@racket[b] is negative --- @racket[+0.0]}
@item{@racket[b] is positive --- @racket[+inf.0]}]}
@item{@racket[(flexpt a b)] where @racket[a] is
negative and @racket[b] is not an integer --- @racket[+nan.0]}
]}
@defproc[(->fl [a exact-integer?]) flonum?]{ @defproc[(->fl [a exact-integer?]) flonum?]{

View File

@ -529,11 +529,71 @@ Returns @racket[(integer-sqrt n)] and @racket[(- n (expt
@defproc[(expt [z number?] [w number?]) number?]{ @defproc[(expt [z number?] [w number?]) number?]{
Returns @racket[z] raised to the power of @racket[w]. If @racket[w] is Returns @racket[z] raised to the power of @racket[w].
exact @racket[0], the result is exact @racket[1]. If @racket[z] is
If @racket[w] is
exact @racket[0], the result is exact @racket[1].
If @racket[w] is @racket[0.0] or @racket[-0.0] and @racket[z] is a @tech{real number}, the
result is @racket[1.0] (even if @racket[z] is @racket[+nan.0]).
If @racket[z] is exact @racket[1], the result is exact @racket[1].
If @racket[z] is @racket[1.0] and @racket[w] is a @tech{real number}, the
result is @racket[1.0] (even if @racket[w] is @racket[+nan.0]).
If @racket[z] is
exact @racket[0] and @racket[w] is negative, the exact @racket[0] and @racket[w] is negative, the
@exnraise[exn:fail:contract:divide-by-zero]. @exnraise[exn:fail:contract:divide-by-zero].
Further special cases when @racket[w] is a @tech{real number}:
@margin-note*{These special cases correspond to @tt{pow} in C99 @cite["C99"],
except when @racket[z] is negative and @racket[w] is a not an
integer.}
@;
@itemlist[#:style 'compact
@item{@racket[(expt 0.0 w)]:
@itemlist[#:style 'compact
@item{@racket[w] is negative --- @racket[+inf.0]}
@item{@racket[w] is positive --- @racket[0.0]}]}
@item{@racket[(expt -0.0 w)]:
@itemlist[#:style 'compact
@item{@racket[w] is negative:
@itemlist[#:style 'compact
@item{@racket[w] is an odd integer --- @racket[-inf.0]}
@item{@racket[w] otherwise rational --- @racket[+inf.0]}]}
@item{@racket[w] is positive:
@itemlist[#:style 'compact
@item{@racket[w] is an odd integer --- @racket[-0.0]}
@item{@racket[w] otherwise rational --- @racket[+0.0]}]}]}
@item{@racket[(expt z -inf.0)] for positive @racket[z]:
@itemlist[#:style 'compact
@item{@racket[z] is less than @racket[1.0] --- @racket[+inf.0]}
@item{@racket[z] is greater than @racket[1.0] --- @racket[+0.0]}]}
@item{@racket[(expt z +inf.0)] for positive @racket[z]:
@itemlist[#:style 'compact
@item{@racket[z] is less than @racket[1.0] --- @racket[+0.0]}
@item{@racket[z] is greater than @racket[1.0] --- @racket[+inf.0]}]}
@item{@racket[(expt -inf.0 w)] for integer @racket[w]:
@itemlist[#:style 'compact
@item{@racket[w] is negative:
@itemlist[#:style 'compact
@item{@racket[w] is odd --- @racket[-0.0]}
@item{@racket[w] is even --- @racket[+0.0]}]}
@item{@racket[w] is positive:
@itemlist[#:style 'compact
@item{@racket[w] is odd --- @racket[-inf.0]}
@item{@racket[w] is even --- @racket[+inf.0]}]}]}
@item{@racket[(expt +inf.0 w)]:
@itemlist[#:style 'compact
@item{@racket[w] is negative --- @racket[+0.0]}
@item{@racket[w] is positive --- @racket[+inf.0]}]}
]
@mz-examples[(expt 2 3) (expt 4 0.5) (expt +inf.0 0)]} @mz-examples[(expt 2 3) (expt 4 0.5) (expt +inf.0 0)]}

View File

@ -79,6 +79,11 @@ The @racketmodname[racket] library combines
@(bibliography @(bibliography
(bib-entry #:key "C99"
#:author "ISO/IEC"
#:title "ISO/IEC 9899:1999 Cor. 3:2007(E))"
#:date "2007")
(bib-entry #:key "Danvy90" (bib-entry #:key "Danvy90"
#:author "Olivier Danvy and Andre Filinski" #:author "Olivier Danvy and Andre Filinski"
#:title "Abstracting Control" #:title "Abstracting Control"

View File

@ -97,6 +97,126 @@
(test-sequence [(2.0 4.0 6.0)] (in-flvector (flvector 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0) 1 6 2)) (test-sequence [(2.0 4.0 6.0)] (in-flvector (flvector 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0) 1 6 2))
(test-sequence [(8.0 6.0 4.0)] (in-flvector (flvector 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0) 7 2 -2)) (test-sequence [(8.0 6.0 4.0)] (in-flvector (flvector 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0) 7 2 -2))
;; ----------------------------------------
;; Check corners of `flexpt':
;; Tests by Neil T.:
(let ()
(define-syntax-rule (check-equal? (flexpt v1 v2) b)
(test b flexpt v1 v2))
;; 2^53 and every larger flonum is even:
(define +big-even.0 (expt 2.0 53))
;; The largest odd flonum:
(define +max-odd.0 (- +big-even.0 1.0))
(define -big-even.0 (- +big-even.0))
(define -max-odd.0 (- +max-odd.0))
(check-equal? (flexpt +0.0 +0.0) +1.0)
(check-equal? (flexpt +0.0 +1.0) +0.0)
(check-equal? (flexpt +0.0 +3.0) +0.0)
(check-equal? (flexpt +0.0 +max-odd.0) +0.0)
(check-equal? (flexpt +0.0 +0.5) +0.0)
(check-equal? (flexpt +0.0 +1.5) +0.0)
(check-equal? (flexpt +0.0 +2.0) +0.0)
(check-equal? (flexpt +0.0 +2.5) +0.0)
(check-equal? (flexpt +0.0 +big-even.0) +0.0)
(check-equal? (flexpt -0.0 +0.0) +1.0)
(check-equal? (flexpt -0.0 +1.0) -0.0)
(check-equal? (flexpt -0.0 +3.0) -0.0)
(check-equal? (flexpt -0.0 +max-odd.0) -0.0)
(check-equal? (flexpt -0.0 +0.5) +0.0)
(check-equal? (flexpt -0.0 +1.5) +0.0)
(check-equal? (flexpt -0.0 +2.0) +0.0)
(check-equal? (flexpt -0.0 +2.5) +0.0)
(check-equal? (flexpt -0.0 +big-even.0) +0.0)
(check-equal? (flexpt +1.0 +0.0) +1.0)
(check-equal? (flexpt +1.0 +0.5) +1.0)
(check-equal? (flexpt +1.0 +inf.0) +1.0)
(check-equal? (flexpt -1.0 +0.0) +1.0)
(check-equal? (flexpt -1.0 +0.5) +nan.0)
(check-equal? (flexpt -1.0 +inf.0) +1.0)
(check-equal? (flexpt +0.5 +inf.0) +0.0)
(check-equal? (flexpt +1.5 +inf.0) +inf.0)
(check-equal? (flexpt +inf.0 +0.0) +1.0)
(check-equal? (flexpt +inf.0 +1.0) +inf.0)
(check-equal? (flexpt +inf.0 +2.0) +inf.0)
(check-equal? (flexpt +inf.0 +inf.0) +inf.0)
(check-equal? (flexpt -inf.0 +0.0) +1.0)
(check-equal? (flexpt -inf.0 +1.0) -inf.0)
(check-equal? (flexpt -inf.0 +3.0) -inf.0)
(check-equal? (flexpt -inf.0 +max-odd.0) -inf.0)
(check-equal? (flexpt -inf.0 +0.5) +inf.0)
(check-equal? (flexpt -inf.0 +1.5) +inf.0)
(check-equal? (flexpt -inf.0 +2.0) +inf.0)
(check-equal? (flexpt -inf.0 +2.5) +inf.0)
(check-equal? (flexpt -inf.0 +big-even.0) +inf.0)
(check-equal? (flexpt -inf.0 +inf.0) +inf.0)
;; Same tests as above, but with negated y
;; This identity should hold for these tests: (flexpt x y) = (/ 1.0 (flexpt x (- y)))
(check-equal? (flexpt +0.0 -0.0) +1.0)
(check-equal? (flexpt +0.0 -1.0) +inf.0)
(check-equal? (flexpt +0.0 -3.0) +inf.0)
(check-equal? (flexpt +0.0 -max-odd.0) +inf.0)
(check-equal? (flexpt +0.0 -0.5) +inf.0)
(check-equal? (flexpt +0.0 -1.5) +inf.0)
(check-equal? (flexpt +0.0 -2.0) +inf.0)
(check-equal? (flexpt +0.0 -2.5) +inf.0)
(check-equal? (flexpt +0.0 -big-even.0) +inf.0)
(check-equal? (flexpt -0.0 -0.0) +1.0)
(check-equal? (flexpt -0.0 -1.0) -inf.0)
(check-equal? (flexpt -0.0 -3.0) -inf.0)
(check-equal? (flexpt -0.0 -max-odd.0) -inf.0)
(check-equal? (flexpt -0.0 -0.5) +inf.0)
(check-equal? (flexpt -0.0 -1.5) +inf.0)
(check-equal? (flexpt -0.0 -2.0) +inf.0)
(check-equal? (flexpt -0.0 -2.5) +inf.0)
(check-equal? (flexpt -0.0 -big-even.0) +inf.0)
(check-equal? (flexpt +1.0 -0.0) +1.0)
(check-equal? (flexpt +1.0 -0.5) +1.0)
(check-equal? (flexpt +1.0 -inf.0) +1.0)
(check-equal? (flexpt -1.0 -0.0) +1.0)
(check-equal? (flexpt -1.0 -0.5) +nan.0)
(check-equal? (flexpt -1.0 -inf.0) +1.0)
(check-equal? (flexpt +0.5 -inf.0) +inf.0)
(check-equal? (flexpt +1.5 -inf.0) +0.0)
(check-equal? (flexpt +inf.0 -0.0) +1.0)
(check-equal? (flexpt +inf.0 -1.0) +0.0)
(check-equal? (flexpt +inf.0 -2.0) +0.0)
(check-equal? (flexpt +inf.0 -inf.0) +0.0)
(check-equal? (flexpt -inf.0 -0.0) +1.0)
(check-equal? (flexpt -inf.0 -1.0) -0.0)
(check-equal? (flexpt -inf.0 -3.0) -0.0)
(check-equal? (flexpt -inf.0 -max-odd.0) -0.0)
(check-equal? (flexpt -inf.0 -0.5) +0.0)
(check-equal? (flexpt -inf.0 -1.5) +0.0)
(check-equal? (flexpt -inf.0 -2.0) +0.0)
(check-equal? (flexpt -inf.0 -2.5) +0.0)
(check-equal? (flexpt -inf.0 -big-even.0) +0.0)
(check-equal? (flexpt -inf.0 -inf.0) +0.0)
;; NaN input
(check-equal? (flexpt +nan.0 +0.0) +1.0)
(check-equal? (flexpt +nan.0 -0.0) +1.0)
(check-equal? (flexpt +1.0 +nan.0) +1.0)
(check-equal? (flexpt -1.0 +nan.0) +nan.0))
;; ----------------------------------------
(report-errs) (report-errs)

View File

@ -1729,18 +1729,14 @@
(test 108.0+29.0i z-round (* 100 (expt 1+i 1/3))) (test 108.0+29.0i z-round (* 100 (expt 1+i 1/3)))
(test 25.0-43.0i z-round (* 100 (expt -8 -1/3))) (test 25.0-43.0i z-round (* 100 (expt -8 -1/3)))
;; This choice doesn't make sense to me, but it fits
;; with other standards and implementations:
(define INF-POWER-OF_NEGATIVE +inf.0)
(test +inf.0 expt 2 +inf.0) (test +inf.0 expt 2 +inf.0)
(test +inf.0 expt +inf.0 10) (test +inf.0 expt +inf.0 10)
(test 0.0 expt +inf.0 -2) (test 0.0 expt +inf.0 -2)
(test 1 expt +inf.0 0) (test 1 expt +inf.0 0)
(test 1.0 expt +inf.0 0.) (test 1.0 expt +inf.0 0.)
(test +inf.0 expt +inf.0 +inf.0) (test +inf.0 expt +inf.0 +inf.0)
(test INF-POWER-OF_NEGATIVE expt -2 +inf.0) (test +nan.0+nan.0i expt -2 +inf.0)
(test INF-POWER-OF_NEGATIVE expt -inf.0 +inf.0) (test +nan.0+nan.0i expt -inf.0 +inf.0)
(test 0.0 expt 2 -inf.0) (test 0.0 expt 2 -inf.0)
(test -inf.0 expt -inf.0 11) (test -inf.0 expt -inf.0 11)
(test +inf.0 expt -inf.0 10) (test +inf.0 expt -inf.0 10)
@ -1749,8 +1745,8 @@
(test 1 expt -inf.0 0) (test 1 expt -inf.0 0)
(test 1.0 expt -inf.0 0.0) (test 1.0 expt -inf.0 0.0)
(test 0.0 expt +inf.0 -inf.0) (test 0.0 expt +inf.0 -inf.0)
(test 0.0 expt -2 -inf.0) (test +nan.0+nan.0i expt -2 -inf.0)
(test 0.0 expt -inf.0 -inf.0) (test +nan.0+nan.0i expt -inf.0 -inf.0)
(test 1 expt +nan.0 0) (test 1 expt +nan.0 0)
(test 0 expt 0 10) (test 0 expt 0 10)
(test 0 expt 0 10.0) (test 0 expt 0 10.0)
@ -1783,9 +1779,9 @@
(test 0 expt 0 1+i) (test 0 expt 0 1+i)
(test 0 expt 0 1-i) (test 0 expt 0 1-i)
(test-nan.0 expt 1.0 +inf.0) (test 1.0 expt 1.0 +inf.0)
(test-nan.0 expt 1.0 -inf.0) (test 1.0 expt 1.0 -inf.0)
(test-nan.0 expt 1.0 +nan.0) (test 1.0 expt 1.0 +nan.0)
(test 0.0 expt 0.0 5) (test 0.0 expt 0.0 5)
(test -0.0 expt -0.0 5) (test -0.0 expt -0.0 5)
@ -1796,13 +1792,13 @@
(test 0.0 expt 0.5 +inf.0) (test 0.0 expt 0.5 +inf.0)
(test +inf.0 expt 0.5 -inf.0) (test +inf.0 expt 0.5 -inf.0)
(test INF-POWER-OF_NEGATIVE expt -0.5 -inf.0) (test +nan.0+nan.0i expt -0.5 -inf.0)
(test +inf.0 expt 1.5 +inf.0) (test +inf.0 expt 1.5 +inf.0)
(test 0.0 expt 1.5 -inf.0) (test 0.0 expt 1.5 -inf.0)
(test 0.0 expt -0.5 +inf.0) (test +nan.0+nan.0i expt -0.5 +inf.0)
(test +inf.0 expt -0.5 -inf.0) (test +nan.0+nan.0i expt -0.5 -inf.0)
(test INF-POWER-OF_NEGATIVE expt -1.5 +inf.0) (test +nan.0+nan.0i expt -1.5 +inf.0)
(test 0.0 expt -1.5 -inf.0) (test +nan.0+nan.0i expt -1.5 -inf.0)
(err/rt-test (expt 0 -1) exn:fail:contract:divide-by-zero?) (err/rt-test (expt 0 -1) exn:fail:contract:divide-by-zero?)
(err/rt-test (expt 0 -1.0) exn:fail:contract:divide-by-zero?) (err/rt-test (expt 0 -1.0) exn:fail:contract:divide-by-zero?)
@ -1821,6 +1817,130 @@
(err/rt-test (expt 1 'a)) (err/rt-test (expt 1 'a))
(err/rt-test (expt 3 'a)) (err/rt-test (expt 3 'a))
;; ----------------------------------------
;; Check corners of `expt':
;; based on the flexpt tests of "flonum.rktl" by Neil T
(let ()
(define-syntax-rule (check-equal? (expt v1 v2) b)
(test b expt v1 v2))
;; 2^53 and every larger flonum is even:
(define +big-even.0 (expt 2.0 53))
;; The largest odd flonum:
(define +max-odd.0 (- +big-even.0 1.0))
(define -big-even.0 (- +big-even.0))
(define -max-odd.0 (- +max-odd.0))
(check-equal? (expt +0.0 +0.0) +1.0)
(check-equal? (expt +0.0 +1.0) +0.0)
(check-equal? (expt +0.0 +3.0) +0.0)
(check-equal? (expt +0.0 +max-odd.0) +0.0)
(check-equal? (expt +0.0 +0.5) +0.0)
(check-equal? (expt +0.0 +1.5) +0.0)
(check-equal? (expt +0.0 +2.0) +0.0)
(check-equal? (expt +0.0 +2.5) +0.0)
(check-equal? (expt +0.0 +big-even.0) +0.0)
(check-equal? (expt -0.0 +0.0) +1.0)
(check-equal? (expt -0.0 +1.0) -0.0)
(check-equal? (expt -0.0 +3.0) -0.0)
(check-equal? (expt -0.0 +max-odd.0) -0.0)
(check-equal? (expt -0.0 +0.5) +0.0)
(check-equal? (expt -0.0 +1.5) +0.0)
(check-equal? (expt -0.0 +2.0) +0.0)
(check-equal? (expt -0.0 +2.5) +0.0)
(check-equal? (expt -0.0 +big-even.0) +0.0)
(check-equal? (expt +1.0 +0.0) +1.0)
(check-equal? (expt +1.0 +0.5) +1.0)
(check-equal? (expt +1.0 +inf.0) +1.0)
(check-equal? (expt -1.0 +0.0) +1.0)
(test 612.0+1e19i 'expt (let ([c (* 1e19 (expt -1.0 +0.5))])
(+ (round (real-part c))
(* 0+1i (round (imag-part c))))))
(check-equal? (expt -1.0 +inf.0) +nan.0+nan.0i)
(check-equal? (expt +0.5 +inf.0) +0.0)
(check-equal? (expt +1.5 +inf.0) +inf.0)
(check-equal? (expt +inf.0 +0.0) +1.0)
(check-equal? (expt +inf.0 +1.0) +inf.0)
(check-equal? (expt +inf.0 +2.0) +inf.0)
(check-equal? (expt +inf.0 +inf.0) +inf.0)
(check-equal? (expt -inf.0 +0.0) +1.0)
(check-equal? (expt -inf.0 +1.0) -inf.0)
(check-equal? (expt -inf.0 +3.0) -inf.0)
(check-equal? (expt -inf.0 +max-odd.0) -inf.0)
(check-equal? (expt -inf.0 +0.5) +inf.0+inf.0i)
(check-equal? (expt -inf.0 +1.5) -inf.0-inf.0i)
(check-equal? (expt -inf.0 +2.0) +inf.0)
(check-equal? (expt -inf.0 +2.5) +inf.0+inf.0i)
(check-equal? (expt -inf.0 +big-even.0) +inf.0)
(check-equal? (expt -inf.0 +inf.0) +nan.0+nan.0i)
;; Same tests as above, but with negated y
;; This identity should hold for these tests: (expt x y) = (/ 1.0 (expt x (- y)))
(check-equal? (expt +0.0 -0.0) +1.0)
(check-equal? (expt +0.0 -1.0) +inf.0)
(check-equal? (expt +0.0 -3.0) +inf.0)
(check-equal? (expt +0.0 -max-odd.0) +inf.0)
(check-equal? (expt +0.0 -0.5) +inf.0)
(check-equal? (expt +0.0 -1.5) +inf.0)
(check-equal? (expt +0.0 -2.0) +inf.0)
(check-equal? (expt +0.0 -2.5) +inf.0)
(check-equal? (expt +0.0 -big-even.0) +inf.0)
(check-equal? (expt -0.0 -0.0) +1.0)
(check-equal? (expt -0.0 -1.0) -inf.0)
(check-equal? (expt -0.0 -3.0) -inf.0)
(check-equal? (expt -0.0 -max-odd.0) -inf.0)
(check-equal? (expt -0.0 -0.5) +inf.0)
(check-equal? (expt -0.0 -1.5) +inf.0)
(check-equal? (expt -0.0 -2.0) +inf.0)
(check-equal? (expt -0.0 -2.5) +inf.0)
(check-equal? (expt -0.0 -big-even.0) +inf.0)
(check-equal? (expt +1.0 -0.0) +1.0)
(check-equal? (expt +1.0 -0.5) +1.0)
(check-equal? (expt +1.0 -inf.0) +1.0)
(check-equal? (expt -1.0 -0.0) +1.0)
(test 612.0-1e19i 'expt (let ([c (* 1e19 (expt -1.0 -0.5))])
(+ (round (real-part c))
(* 0+1i (round (imag-part c))))))
(check-equal? (expt -1.0 -inf.0) +nan.0+nan.0i)
(check-equal? (expt +0.5 -inf.0) +inf.0)
(check-equal? (expt +1.5 -inf.0) +0.0)
(check-equal? (expt +inf.0 -0.0) +1.0)
(check-equal? (expt +inf.0 -1.0) +0.0)
(check-equal? (expt +inf.0 -2.0) +0.0)
(check-equal? (expt +inf.0 -inf.0) +0.0)
(check-equal? (expt -inf.0 -0.0) +1.0)
(check-equal? (expt -inf.0 -1.0) -0.0)
(check-equal? (expt -inf.0 -3.0) -0.0)
(check-equal? (expt -inf.0 -max-odd.0) -0.0)
(check-equal? (expt -inf.0 -0.5) +0.0-0.0i)
(check-equal? (expt -inf.0 -1.5) -0.0+0.0i)
(check-equal? (expt -inf.0 -2.0) +0.0)
(check-equal? (expt -inf.0 -2.5) 0.0-0.0i)
(check-equal? (expt -inf.0 -big-even.0) +0.0)
(check-equal? (expt -inf.0 -inf.0) +nan.0+nan.0i)
;; NaN input
(check-equal? (expt +nan.0 +0.0) +1.0)
(check-equal? (expt +nan.0 -0.0) +1.0)
(check-equal? (expt +1.0 +nan.0) +1.0)
(check-equal? (expt -1.0 +nan.0) +nan.0+nan.0i))
;;;;From: fred@sce.carleton.ca (Fred J Kaudel) ;;;;From: fred@sce.carleton.ca (Fred J Kaudel)
;;; Modified by jaffer. ;;; Modified by jaffer.
(define f3.9 (string->number "3.9")) (define f3.9 (string->number "3.9"))

View File

@ -408,8 +408,8 @@
(bin-exact -0.125 'flexpt -2.0 -3.0 #t) (bin-exact -0.125 'flexpt -2.0 -3.0 #t)
(bin-exact +nan.0 'flexpt -1.0 3.1 #t) (bin-exact +nan.0 'flexpt -1.0 3.1 #t)
(bin-exact 0.0 'flexpt 0.0 10.0 #t) (bin-exact 0.0 'flexpt 0.0 10.0 #t)
(bin-exact +nan.0 'flexpt 0.0 -1.0 #t) (bin-exact +inf.0 'flexpt 0.0 -1.0 #t)
(bin-exact +nan.0 'flexpt 0.0 0.0 #t) (bin-exact +1.0 'flexpt 0.0 0.0 #t)
(bin-exact +nan.0 'flexpt +nan.0 2.7 #t) (bin-exact +nan.0 'flexpt +nan.0 2.7 #t)
(bin-exact +nan.0 'flexpt 2.7 +nan.0 #t) (bin-exact +nan.0 'flexpt 2.7 +nan.0 #t)
(bin-exact +nan.0 'flexpt +nan.0 +nan.0 #t) (bin-exact +nan.0 'flexpt +nan.0 +nan.0 #t)

View File

@ -82,7 +82,6 @@
# include "uconfig.h" # include "uconfig.h"
# define USE_EXPLICT_FP_FORM_CHECK # define USE_EXPLICT_FP_FORM_CHECK
# define POW_HANDLES_INF_CORRECTLY
# include <errno.h> # include <errno.h>
# ifdef ECHRNG # ifdef ECHRNG
@ -144,8 +143,6 @@
# define SELECT_INCLUDE # define SELECT_INCLUDE
# define POW_HANDLES_INF_CORRECTLY
# define USE_FCNTL_O_NONBLOCK # define USE_FCNTL_O_NONBLOCK
# define USE_TIMEZONE_VAR_W_DLS # define USE_TIMEZONE_VAR_W_DLS
@ -261,7 +258,6 @@
#endif #endif
# define USE_IEEE_FP_PREDS # define USE_IEEE_FP_PREDS
# define POW_HANDLES_INF_CORRECTLY
# define USE_DYNAMIC_FDSET_SIZE # define USE_DYNAMIC_FDSET_SIZE
@ -311,7 +307,6 @@
#endif #endif
# define USE_IEEE_FP_PREDS # define USE_IEEE_FP_PREDS
# define POW_HANDLES_INF_CORRECTLY
# define USE_DYNAMIC_FDSET_SIZE # define USE_DYNAMIC_FDSET_SIZE
@ -382,9 +377,6 @@
# define USE_UNDERSCORE_SETJMP # define USE_UNDERSCORE_SETJMP
# define USE_IEEE_FP_PREDS # define USE_IEEE_FP_PREDS
# ifndef ASM_DBLPREC_CONTROL_87
# define POW_HANDLES_INF_CORRECTLY
# endif
# define USE_DYNAMIC_FDSET_SIZE # define USE_DYNAMIC_FDSET_SIZE
@ -773,6 +765,7 @@
#ifdef MZ_NO_JIT_SSE #ifdef MZ_NO_JIT_SSE
# define ASM_DBLPREC_CONTROL_87 # define ASM_DBLPREC_CONTROL_87
#endif #endif
# define POW_HANDLES_CASES_CORRECTLY
# define MZ_JIT_USE_MPROTECT # define MZ_JIT_USE_MPROTECT
@ -861,7 +854,6 @@
# define NO_STAT_PROC # define NO_STAT_PROC
# define DONT_IGNORE_PIPE_SIGNAL # define DONT_IGNORE_PIPE_SIGNAL
# define POW_HANDLES_INF_CORRECTLY
# define TRIG_ZERO_NEEDS_SIGN_CHECK # define TRIG_ZERO_NEEDS_SIGN_CHECK
# define MACOS_UNICODE_SUPPORT # define MACOS_UNICODE_SUPPORT
@ -947,7 +939,6 @@
# define DONT_IGNORE_PIPE_SIGNAL # define DONT_IGNORE_PIPE_SIGNAL
# define DONT_IGNORE_FPE_SIGNAL # define DONT_IGNORE_FPE_SIGNAL
# define POW_HANDLES_INF_CORRECTLY
# define USE_PALM_INF_TESTS # define USE_PALM_INF_TESTS
# define FLAGS_ALREADY_SET # define FLAGS_ALREADY_SET
@ -1323,9 +1314,9 @@
SunOS/Solaris, and HP/UX by explicit pre-checking the form of the SunOS/Solaris, and HP/UX by explicit pre-checking the form of the
number and looking for values that are obviously +inf.0 or -inf.0 */ number and looking for values that are obviously +inf.0 or -inf.0 */
/* POW_HANDLES_INF_CORRECTLY indicates that thw pow() library procedure /* POW_HANDLES_CASES_CORRECTLY indicates that the pow() library procedure
handles +/-inf.0 correctly. Otherwise, code in inserted to specifically handles all +/-inf.0, +/-0.0, or +nan.0 cases according to C99. This
check for infinite arguments. */ might save time on redundant checking before Racket calls pow(). */
/* ATAN2_DOESNT_WORK_WITH_INFINITIES indicates that atan2(+/-inf, +/-inf) /* ATAN2_DOESNT_WORK_WITH_INFINITIES indicates that atan2(+/-inf, +/-inf)
is not the same as atan2(1, 1). */ is not the same as atan2(1, 1). */

View File

@ -300,6 +300,7 @@ Scheme_Object *scheme_complex_power(const Scheme_Object *base, const Scheme_Obje
Scheme_Complex *cb = (Scheme_Complex *)base; Scheme_Complex *cb = (Scheme_Complex *)base;
Scheme_Complex *ce = (Scheme_Complex *)exponent; Scheme_Complex *ce = (Scheme_Complex *)exponent;
double a, b, c, d, bm, ba, nm, na, r1, r2; double a, b, c, d, bm, ba, nm, na, r1, r2;
int d_is_zero;
if ((ce->i == zero) && !SCHEME_FLOATP(ce->r)) { if ((ce->i == zero) && !SCHEME_FLOATP(ce->r)) {
if (SCHEME_INTP(ce->r) || SCHEME_BIGNUMP(ce->r)) if (SCHEME_INTP(ce->r) || SCHEME_BIGNUMP(ce->r))
@ -310,13 +311,17 @@ Scheme_Object *scheme_complex_power(const Scheme_Object *base, const Scheme_Obje
b = scheme_get_val_as_double(cb->i); b = scheme_get_val_as_double(cb->i);
c = scheme_get_val_as_double(ce->r); c = scheme_get_val_as_double(ce->r);
d = scheme_get_val_as_double(ce->i); d = scheme_get_val_as_double(ce->i);
d_is_zero = (ce->i == zero);
bm = sqrt(a * a + b * b); bm = sqrt(a * a + b * b);
ba = atan2(b, a); ba = atan2(b, a);
/* New mag & angle */ /* New mag & angle */
nm = pow(bm, c) * exp(-(ba * d)); nm = scheme_double_expt(bm, c) * exp(-(ba * d));
na = log(bm) * d + ba * c; if (d_is_zero) /* precision here can avoid NaNs */
na = ba * c;
else
na = log(bm) * d + ba * c;
r1 = nm * cos(na); r1 = nm * cos(na);
r2 = nm * sin(na); r2 = nm * sin(na);

View File

@ -455,9 +455,9 @@ static Scheme_Object *unary_minus(const Scheme_Object *n)
#define ret_1other(n1, n2) if (SAME_OBJ(n1, scheme_make_integer(1))) return (Scheme_Object *)n2 #define ret_1other(n1, n2) if (SAME_OBJ(n1, scheme_make_integer(1))) return (Scheme_Object *)n2
#define ret_zero(n1, n2) if (SAME_OBJ(n1, scheme_make_integer(0))) return scheme_make_integer(0) #define ret_zero(n1, n2) if (SAME_OBJ(n1, scheme_make_integer(0))) return scheme_make_integer(0)
GEN_BIN_OP(scheme_bin_plus, "+", ADD, F_ADD, FS_ADD, scheme_bignum_add, scheme_rational_add, scheme_complex_add, GEN_RETURN_N2, GEN_RETURN_N1, NO_NAN_CHECK, NO_NAN_CHECK, ret_other, cx_NO_CHECK, ret_other, cx_NO_CHECK) GEN_BIN_OP(scheme_bin_plus, "+", ADD, F_ADD, FS_ADD, scheme_bignum_add, scheme_rational_add, scheme_complex_add, GEN_RETURN_N2, GEN_RETURN_N1, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, ret_other, cx_NO_CHECK, ret_other, cx_NO_CHECK)
GEN_BIN_OP(scheme_bin_minus, "-", SUBTRACT, F_SUBTRACT, FS_SUBTRACT, scheme_bignum_subtract, scheme_rational_subtract, scheme_complex_subtract, GEN_SINGLE_SUBTRACT_N2, GEN_RETURN_N1, NO_NAN_CHECK, NO_NAN_CHECK, cx_NO_CHECK, cx_NO_CHECK, ret_other, cx_NO_CHECK) GEN_BIN_OP(scheme_bin_minus, "-", SUBTRACT, F_SUBTRACT, FS_SUBTRACT, scheme_bignum_subtract, scheme_rational_subtract, scheme_complex_subtract, GEN_SINGLE_SUBTRACT_N2, GEN_RETURN_N1, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, cx_NO_CHECK, cx_NO_CHECK, ret_other, cx_NO_CHECK)
GEN_BIN_OP(scheme_bin_mult, "*", MULTIPLY, F_MULTIPLY, FS_MULTIPLY, scheme_bignum_multiply, scheme_rational_multiply, scheme_complex_multiply, GEN_RETURN_0, GEN_RETURN_0, NO_NAN_CHECK, NO_NAN_CHECK, ret_zero, ret_1other, ret_zero, ret_1other) GEN_BIN_OP(scheme_bin_mult, "*", MULTIPLY, F_MULTIPLY, FS_MULTIPLY, scheme_bignum_multiply, scheme_rational_multiply, scheme_complex_multiply, GEN_RETURN_0, GEN_RETURN_0, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, NO_NAN_CHECK, ret_zero, ret_1other, ret_zero, ret_1other)
GEN_BIN_DIV_OP(scheme_bin_div, "/", DIVIDE, F_DIVIDE, FS_DIVIDE, scheme_make_rational, scheme_rational_divide, scheme_complex_divide, ret_zero, cx_NO_CHECK, cx_NO_CHECK, ret_1other) GEN_BIN_DIV_OP(scheme_bin_div, "/", DIVIDE, F_DIVIDE, FS_DIVIDE, scheme_make_rational, scheme_rational_divide, scheme_complex_divide, ret_zero, cx_NO_CHECK, cx_NO_CHECK, ret_1other)
GEN_NARY_OP(static, plus, "+", scheme_bin_plus, 0, SCHEME_NUMBERP, "number?", GEN_IDENT) GEN_NARY_OP(static, plus, "+", scheme_bin_plus, 0, SCHEME_NUMBERP, "number?", GEN_IDENT)

View File

@ -1278,23 +1278,29 @@ rational_p(int argc, Scheme_Object *argv[])
return (is_rational(argv[0]) ? scheme_true : scheme_false); return (is_rational(argv[0]) ? scheme_true : scheme_false);
} }
XFORM_NONGCING static int double_is_integer(double d)
{
# ifdef NAN_EQUALS_ANYTHING
if (MZ_IS_NAN(d))
return 0;
# endif
if (MZ_IS_INFINITY(d))
return 0;
if (floor(d) == d)
return 1;
return 0;
}
int scheme_is_integer(const Scheme_Object *o) int scheme_is_integer(const Scheme_Object *o)
{ {
if (SCHEME_INTP(o) || SCHEME_BIGNUMP(o)) if (SCHEME_INTP(o) || SCHEME_BIGNUMP(o))
return 1; return 1;
if (SCHEME_FLOATP(o)) { if (SCHEME_FLOATP(o))
double d; return double_is_integer(SCHEME_FLOAT_VAL(o));
d = SCHEME_FLOAT_VAL(o);
# ifdef NAN_EQUALS_ANYTHING
if (MZ_IS_NAN(d))
return 0;
# endif
if (MZ_IS_INFINITY(d))
return 0;
if (floor(d) == d)
return 1;
}
return 0; return 0;
} }
@ -2509,82 +2515,111 @@ static Scheme_Object *fixnum_expt(intptr_t x, intptr_t y)
} }
} }
#ifdef POW_HANDLES_INF_CORRECTLY #ifdef ASM_DBLPREC_CONTROL_87
# define sch_pow pow static double protected_pow(double x, double y)
{
/* libm's pow() implementation seems to sometimes rely on
extended precision in pow(), so reset the control
word while calling pow(); note that the x87 control
word is thread-specific */
to_extended_prec();
x = pow(x, y);
to_double_prec();
return x;
}
#else
# define protected_pow pow
#endif
#ifdef POW_HANDLES_CASES_CORRECTLY
# define sch_pow protected_pow
#else #else
static double sch_pow(double x, double y) static double sch_pow(double x, double y)
{ {
if (MZ_IS_POS_INFINITY(y)) { /* Explciitly handle all cases described by C99 */
if ((x == 1.0) || (x == -1.0)) if (x == 1.0)
return not_a_number_val; return 1.0; /* even for NaN */
else if (y == 0.0)
return 1.0; /* even for NaN */
else if (MZ_IS_NAN(x))
return not_a_number_val;
else if (MZ_IS_NAN(y))
return not_a_number_val;
else if (x == 0.0) {
int neg = 0;
if (y < 0) {
neg = 1;
y = -y;
}
if (fmod(y, 2.0) == 1.0) {
if (neg) {
if (minus_zero_p(x))
return scheme_minus_infinity_val;
else
return scheme_infinity_val;
} else
return x;
} else {
if (neg)
return scheme_infinity_val;
else
return 0.0;
}
} else if (MZ_IS_POS_INFINITY(y)) {
if (x == -1.0)
return 1.0;
else if ((x < 1.0) && (x > -1.0)) else if ((x < 1.0) && (x > -1.0))
return 0.0; return 0.0;
else else
return scheme_infinity_val; return scheme_infinity_val;
} else if (MZ_IS_NEG_INFINITY(y)) { } else if (MZ_IS_NEG_INFINITY(y)) {
if ((x == 1.0) || (x == -1.0)) if (x == -1.0)
return not_a_number_val; return 1.0;
else if ((x < 1.0) && (x > -1.0)) else if ((x < 1.0) && (x > -1.0))
return scheme_infinity_val; return scheme_infinity_val;
else else
return 0.0; return 0.0;
} else if (MZ_IS_POS_INFINITY(x)) { } else if (MZ_IS_POS_INFINITY(x)) {
if (y == 0.0) if (y < 0)
return 1.0;
else if (y < 0)
return 0.0; return 0.0;
else else
return scheme_infinity_val; return scheme_infinity_val;
} else if (MZ_IS_NEG_INFINITY(x)) { } else if (MZ_IS_NEG_INFINITY(x)) {
if (y == 0.0) int neg = 0;
return 1.0; if (y < 0) {
else { neg = 1;
int neg = 0; y = -y;
if (y < 0) { }
neg = 1; if (fmod(y, 2.0) == 1.0) {
y = -y; if (neg)
} return scheme_floating_point_nzero;
if (fmod(y, 2.0) == 1.0) { else
if (neg) return scheme_minus_infinity_val;
return scheme_floating_point_nzero; } else {
else if (neg)
return scheme_minus_infinity_val; return 0.0;
} else { else
if (neg) return scheme_infinity_val;
return 0.0;
else
return scheme_infinity_val;
}
} }
} else { } else {
#ifdef ASM_DBLPREC_CONTROL_87 return protected_pow(x, y);
/* libm's pow() implementation seems to rely on
extended precision in pow(), so reset the control
word while calling pow(); note that the x87 control
word is thread-specific */
to_extended_prec();
#endif
x = pow(x, y);
#ifdef ASM_DBLPREC_CONTROL_87
to_double_prec();
#endif
return x;
} }
} }
#endif #endif
GEN_BIN_PROT(bin_expt); GEN_BIN_PROT(bin_expt);
# define F_EXPT(x, y) (((x < 0.0) && (y != floor(y))) \ # define F_EXPT(x, y) (((x < 0.0) && !double_is_integer(y)) \
? scheme_complex_power(scheme_real_to_complex(scheme_make_double(x)), \ ? scheme_complex_power(scheme_real_to_complex(scheme_make_double(x)), \
scheme_real_to_complex(scheme_make_double(y))) \ scheme_real_to_complex(scheme_make_double(y))) \
: scheme_make_double(sch_pow((double)x, (double)y))) : scheme_make_double(sch_pow((double)x, (double)y)))
# define FS_EXPT(x, y) (((x < 0.0) && (y != floor(y))) \ # define FS_EXPT(x, y) (((x < 0.0) && !double_is_integer(y)) \
? scheme_complex_power(scheme_real_to_complex(scheme_make_float(x)), \ ? scheme_complex_power(scheme_real_to_complex(scheme_make_float(x)), \
scheme_real_to_complex(scheme_make_float(y))) \ scheme_real_to_complex(scheme_make_float(y))) \
: scheme_make_float(sch_pow((double)x, (double)y))) : scheme_make_float(sch_pow((double)x, (double)y)))
static GEN_BIN_OP(bin_expt, "expt", fixnum_expt, F_EXPT, FS_EXPT, scheme_generic_integer_power, scheme_rational_power, scheme_complex_power, GEN_RETURN_0_USUALLY, GEN_RETURN_1, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK) static GEN_BIN_OP(bin_expt, "expt", fixnum_expt, F_EXPT, FS_EXPT, scheme_generic_integer_power, scheme_rational_power, scheme_complex_power, GEN_RETURN_0_USUALLY, GEN_RETURN_1, NO_NAN_CHECK, NO_NAN_CHECK, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK)
Scheme_Object * Scheme_Object *
scheme_expt(int argc, Scheme_Object *argv[]) scheme_expt(int argc, Scheme_Object *argv[])
@ -2731,14 +2766,7 @@ scheme_expt(int argc, Scheme_Object *argv[])
} }
double scheme_double_expt(double x, double y) { double scheme_double_expt(double x, double y) {
if ((x < 0) && (floor(y) != y)) return sch_pow(x, y);
return not_a_number_val;
else if ((x == 0.0) && (y <= 0))
return not_a_number_val;
else if (MZ_IS_NAN(x) || MZ_IS_NAN(y))
return not_a_number_val;
else
return sch_pow(x, y);
} }
Scheme_Object *scheme_checked_make_rectangular (int argc, Scheme_Object *argv[]) Scheme_Object *scheme_checked_make_rectangular (int argc, Scheme_Object *argv[])

View File

@ -506,8 +506,8 @@ negative_p (int argc, Scheme_Object *argv[])
#define MAX_IZI(a, b) bin_max(IZI_REAL_PART(a), IZI_REAL_PART(b)) #define MAX_IZI(a, b) bin_max(IZI_REAL_PART(a), IZI_REAL_PART(b))
#define MIN_IZI(a, b) bin_min(IZI_REAL_PART(a), IZI_REAL_PART(b)) #define MIN_IZI(a, b) bin_min(IZI_REAL_PART(a), IZI_REAL_PART(b))
static GEN_BIN_OP(bin_max, "max", MAX, F_MAX, FS_MAX, scheme_bignum_max, scheme_rational_max, MAX_IZI, GEN_OMIT, GEN_OMIT, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK) static GEN_BIN_OP(bin_max, "max", MAX, F_MAX, FS_MAX, scheme_bignum_max, scheme_rational_max, MAX_IZI, GEN_OMIT, GEN_OMIT, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK)
static GEN_BIN_OP(bin_min, "min", MIN, F_MIN, FS_MIN, scheme_bignum_min, scheme_rational_min, MIN_IZI, GEN_OMIT, GEN_OMIT, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK) static GEN_BIN_OP(bin_min, "min", MIN, F_MIN, FS_MIN, scheme_bignum_min, scheme_rational_min, MIN_IZI, GEN_OMIT, GEN_OMIT, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, NAN_RETURNS_NAN, NAN_RETURNS_SNAN, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK, cx_NO_CHECK)
GEN_TWOARY_OP(static, sch_max, "max", bin_max, SCHEME_REALP, "real?") GEN_TWOARY_OP(static, sch_max, "max", bin_max, SCHEME_REALP, "real?")
GEN_TWOARY_OP(static, sch_min, "min", bin_min, SCHEME_REALP, "real?") GEN_TWOARY_OP(static, sch_min, "min", bin_min, SCHEME_REALP, "real?")

View File

@ -540,14 +540,14 @@ name (const Scheme_Object *n1, const Scheme_Object *n2) \
# define NAN_CHECK_0(x) if (MZ_IS_NAN(x)) return 0 # define NAN_CHECK_0(x) if (MZ_IS_NAN(x)) return 0
#define GEN_BIN_OP(name, scheme_name, iop, fop, fsop, bn_op, rop, cxop, exzeopl, exzeopr, nanckop, snanckop, c0_1, c1_1, c0_2, c1_2) \ #define GEN_BIN_OP(name, scheme_name, iop, fop, fsop, bn_op, rop, cxop, exzeopl, exzeopr, nanckop, snanckop, nanckop_more, snanckop_more, c0_1, c1_1, c0_2, c1_2) \
GEN_BIN_THING(Scheme_Object *, name, scheme_name, \ GEN_BIN_THING(Scheme_Object *, name, scheme_name, \
iop, fop, fsop, bn_op, rop, cxop, \ iop, fop, fsop, bn_op, rop, cxop, \
GEN_OMIT, GEN_FIRST_ONLY, \ GEN_OMIT, GEN_FIRST_ONLY, \
0, 0, 0, 0, \ 0, 0, 0, 0, \
0, 0, 0, 0, \ 0, 0, 0, 0, \
GEN_SCHEME_BOOL_APPLY, badfunc, badfunc, badfunc, badfunc, \ GEN_SCHEME_BOOL_APPLY, badfunc, badfunc, badfunc, badfunc, \
nanckop, snanckop, nanckop, snanckop, \ nanckop, snanckop, nanckop_more, snanckop_more, \
GEN_IDENT, GEN_IDENT, exzeopl, exzeopr, "number?", GEN_TOI, \ GEN_IDENT, GEN_IDENT, exzeopl, exzeopr, "number?", GEN_TOI, \
c0_1, c1_1, c0_2, c1_2) c0_1, c1_1, c0_2, c1_2)