fix broken exponential asm

infinities were getting converted into nans. the new code simply tests
for infinity and replaces it with a large magnitude value of the same
sign.

also, the fcomi instruction is apparently not part of the i387
instruction set, so avoid using it.

diff --git a/src/math/i386/exp.s b/src/math/i386/exp.s
index 18f6cd67f6fbc5f24056ba6e9350ffdaec666092..f4769d59c1ee499beacc8df5974ab45be60b9c8b 100644 (file)
--- a/src/math/i386/exp.s
+++ b/src/math/i386/exp.s
@@ -34,6 +34,15 @@ exp:
 .type exp2,@function
 exp2:
        fldl 4(%esp)
+1:     fxam
+       fnstsw %ax
+       sahf
+       jnp 1f
+       jnc 1f
+       fstps 4(%esp)
+       mov $0xfe,%al
+       and %al,7(%esp)
+       flds 4(%esp)
 1:     fld %st(0)
        frndint
        fxch %st(1)

diff --git a/src/math/i386/expm1.s b/src/math/i386/expm1.s
index d6d511ef4689e3f0271d1cf67fa5dca6167f98d3..bbb5d12eae7dee1715751e64fdf82c0e5bd1ff8f 100644 (file)
--- a/src/math/i386/expm1.s
+++ b/src/math/i386/expm1.s
@@ -14,12 +14,23 @@ expm1l:
 .type expm1,@function
 expm1:
        fldl 4(%esp)
+1:     fxam
+       fnstsw %ax
+       sahf
+       jnp 1f
+       jnc 1f
+       fstps 4(%esp)
+       mov $0xfe,%al
+       and %al,7(%esp)
+       flds 4(%esp)
 1:     fldl2e
        fmulp
        fld %st(0)
        frndint
        fldz
-       fcompi
+       fcomp
+       fnstsw %ax
+       sahf
        jnz 1f
        fstp %st(0)
        f2xm1