Joseph Ferano
aeeb6de59d
The prelude's declare surface was five f32 functions, and the five were there because somebody needed each one. Everything else a caller wanted was written as a declare at the top of their own file -- the identical libm call with none of the caveats written down. So the rest of libm is here: tan, the three inverses, the three logarithms, exp, fmod, hypot, cbrt, fabs, and an f64 face for every one of them including the five that already existed. A declare is a line, a symbol already on the link, and nothing in either backend, which is why this was cheap enough to do completely rather than one function at a time. The f64 half is not decoration. f32 is what a position is; f64 is what a measurement is -- the clock, parse-f64, format-f64, any sum over more than a few thousand terms -- and having only the f32 face forced a cast down and back at each of those boundaries, which is where the precision went. The paragraph the sqrt note draws for itself is now drawn once for the family: IEEE-754 specifies sqrt, fabs, floor, ceil, round and fmod as exact or correctly rounded, so those agree bit for bit across glibc, musl and wasi-libc; it requires nothing of the rest, so the sand-grid rule covers all of them unchanged. floor, ceil and round are Flan at f32 and libm at f64, and that is not an inconsistency: the f32 bodies work because every f32 with a fraction fits in an i32, and at f64 that trick is gone. abs-i32 and abs-i64 are Flan, one per width because min and max are builtins and no generic covers the numeric types. pi and tau at both widths, written out rather than derived so the compiler rounds each literal once. programs/math3.flan covers it at values that are exact in binary, so nothing pins one libm's last bit. The -O0 case is the one that matters: at -O2 LLVM folds a call over two literals and leaves no symbol to resolve, which is how a missing -lm hid the first time.
93 lines
3.8 KiB
Plaintext
93 lines
3.8 KiB
Plaintext
;;;; The rest of libm, at both widths — what math.flan and math2.flan left out.
|
|
;;;;
|
|
;;;; The rule those two set holds here unchanged and is the only reason this
|
|
;;;; file looks the way it does: none of these is correctly rounded under
|
|
;;;; IEEE-754 except sqrt, fabs, the rounding three and fmod, so every value
|
|
;;;; below is one whose answer is exact in binary — zero, one, a power of two,
|
|
;;;; a perfect square, a perfect cube. A case that pinned glibc's last bit
|
|
;;;; would pass here and fail on wasi-libc.
|
|
;;;;
|
|
;;;; The acceptance table builds this at -O0 as well, and that run is the one
|
|
;;;; that matters: at -O2 LLVM folds a libm call over two literals and leaves
|
|
;;;; no symbol to resolve, which is how a missing -lm hid the first time.
|
|
|
|
(defn show [x f32] ()
|
|
(print x)
|
|
(print " "))
|
|
|
|
(defn show64 [x f64] ()
|
|
(print x)
|
|
(print " "))
|
|
|
|
(defn main [] i32
|
|
;; The f32 half. tan, the three inverses, the three logarithms and exp.
|
|
(show (tan-f32 0.0)) ; 0
|
|
(show (asin-f32 0.0)) ; 0
|
|
(show (acos-f32 1.0)) ; 0
|
|
(show (atan-f32 0.0)) ; 0
|
|
(show (log-f32 1.0)) ; 0
|
|
(show (log2-f32 8.0)) ; 3
|
|
(show (log10-f32 1000.0)) ; 3
|
|
(show (exp-f32 0.0)) ; 1
|
|
(println "")
|
|
|
|
(show (fmod-f32 7.0 4.0)) ; 3
|
|
;; The sign follows the dividend and not the divisor, which is the line a
|
|
;; caller reaching for a modulo gets wrong. Written down as a case.
|
|
(show (fmod-f32 -1.0 3.0)) ; -1
|
|
(show (hypot-f32 3.0 4.0)) ; 5
|
|
(show (cbrt-f32 27.0)) ; 3
|
|
;; The negative root, where (pow-f32 x 0.33333334) would be NaN: pow goes
|
|
;; through a logarithm and cbrt does not.
|
|
(show (cbrt-f32 -8.0)) ; -2
|
|
(show (abs-f32 -2.5)) ; 2.5
|
|
(println "")
|
|
|
|
;; The f64 half, at the same exact values. This block is the whole point of
|
|
;; the f64 face existing: before it, every one of these was a cast down to
|
|
;; f32 and back, and the cast down is where the precision went.
|
|
(show64 (sqrt-f64 16.0)) ; 4
|
|
(show64 (sin-f64 0.0)) ; 0
|
|
(show64 (cos-f64 0.0)) ; 1
|
|
(show64 (tan-f64 0.0)) ; 0
|
|
(show64 (asin-f64 0.0)) ; 0
|
|
(show64 (acos-f64 1.0)) ; 0
|
|
(show64 (atan-f64 0.0)) ; 0
|
|
(show64 (atan2-f64 0.0 1.0)) ; 0
|
|
(println "")
|
|
|
|
(show64 (log-f64 1.0)) ; 0
|
|
(show64 (log2-f64 1024.0)) ; 10
|
|
(show64 (log10-f64 100.0)) ; 2
|
|
(show64 (exp-f64 0.0)) ; 1
|
|
(show64 (pow-f64 2.0 10.0)) ; 1024
|
|
(show64 (fmod-f64 7.0 4.0)) ; 3
|
|
(show64 (hypot-f64 3.0 4.0)) ; 5
|
|
(show64 (cbrt-f64 8.0)) ; 2
|
|
(show64 (abs-f64 -1.5)) ; 1.5
|
|
(println "")
|
|
|
|
;; The f64 rounding family, which is libm's where the f32 one is Flan's —
|
|
;; and it agrees with the Flan one where they overlap: half away from zero,
|
|
;; so -2.5 goes to -3 and not to -2.
|
|
(show64 (floor-f64 -2.5)) ; -3
|
|
(show64 (ceil-f64 -2.5)) ; -2
|
|
(show64 (round-f64 -2.5)) ; -3
|
|
(show64 (round-f64 2.5)) ; 3
|
|
(println "")
|
|
|
|
;; Integer magnitude, one per width.
|
|
(print (abs-i32 -7)) (print " ") ; 7
|
|
(print (abs-i64 (i64 -7))) (print " ") ; 7
|
|
(print (abs-i32 7)) (print " ") ; 7
|
|
;; tau is 2pi at both widths. Pinning the relation rather than the digits is
|
|
;; what catches a constant written to too few of them.
|
|
(print (= tau-f32 (* 2.0 pi-f32))) (print " ")
|
|
(print (= tau-f64 (* 2.0 pi-f64)))
|
|
(println "")
|
|
|
|
;; pi is the one value here that can be pinned without pinning a libm: it is
|
|
;; a literal the compiler rounds, so it is the same on every target.
|
|
(println (and (> pi-f64 3.14159265) (< pi-f64 3.14159266)))
|
|
0)
|