diff options
| author | nsz <nsz@port70.net> | 2012-08-08 20:18:16 +0200 |
|---|---|---|
| committer | nsz <nsz@port70.net> | 2012-08-08 20:18:16 +0200 |
| commit | 1fb01691547e85be253a4f1543be3920de651e27 (patch) | |
| tree | c7957da95fa496f77fea7e118faa1440003825ce /src/math | |
| parent | ae0b9da48c91087c5ab78e4918deb69665d0ccc6 (diff) | |
| download | musl-1fb01691547e85be253a4f1543be3920de651e27.tar.gz musl-1fb01691547e85be253a4f1543be3920de651e27.tar.xz musl-1fb01691547e85be253a4f1543be3920de651e27.zip | |
math: fix exp.s on i386 and x86_64 so the exception flags are correct
exp(inf), exp(-inf), exp(nan) used to raise wrong flags
Diffstat (limited to 'src/math')
| -rw-r--r-- | src/math/i386/exp.s | 39 | ||||
| -rw-r--r-- | src/math/x86_64/expl.s | 36 |
2 files changed, 35 insertions, 40 deletions
diff --git a/src/math/i386/exp.s b/src/math/i386/exp.s index ca0de1d4..c7f5ad0f 100644 --- a/src/math/i386/exp.s +++ b/src/math/i386/exp.s @@ -69,18 +69,18 @@ exp: exp2: fldl 4(%esp) 1: pushl $0x467ff000 - flds (%esp) + flds (%esp) # 16380 xorl %eax,%eax pushl $0x80000000 push %eax fld %st(1) fabs - fucom %st(1) + fucomp %st(1) fnstsw fstp %st(0) - fstp %st(0) sahf - ja 2f + ja 3f # |x| > 16380 + jp 2f # x is nan (avoid invalid except in fistp) fld %st(0) fistpl 8(%esp) fildl 8(%esp) @@ -90,30 +90,27 @@ exp2: add %eax,8(%esp) f2xm1 fld1 - faddp - fldt (%esp) + faddp # 2^(x-rint(x)) + fldt (%esp) # 2^rint(x) fmulp fstp %st(1) - add $12,%esp +2: add $12,%esp ret -2: fld %st(0) +3: fld %st(0) fstpt (%esp) - mov 9(%esp),%ah - and $0x7f,%ah - cmp $0x7f,%ah - jne 1f - decb 9(%esp) - fstp %st(0) - fldt (%esp) -1: fld %st(0) + fld1 + mov 8(%esp),%ax + and $0x7fff,%ax + cmp $0x7fff,%ax + je 1f # x = +-inf + fld %st(1) frndint - fxch %st(1) - fsub %st(1) + fxch %st(2) + fsub %st(2) # st(0)=x-rint(x), st(1)=1, st(2)=rint(x) f2xm1 - fld1 - faddp - fscale + faddp # 2^(x-rint(x)) +1: fscale fstp %st(1) add $12,%esp ret diff --git a/src/math/x86_64/expl.s b/src/math/x86_64/expl.s index 64c1c78e..740bc779 100644 --- a/src/math/x86_64/expl.s +++ b/src/math/x86_64/expl.s @@ -40,7 +40,7 @@ exp2l: mov %eax,-20(%rsp) xor %eax,%eax mov %eax,-24(%rsp) - flds -16(%rsp) + flds -16(%rsp) # 16380 fld %st(1) fabs fucom %st(1) @@ -48,7 +48,8 @@ exp2l: fstp %st(0) fstp %st(0) sahf - ja 2f + ja 3f # |x| > 16380 + jp 2f # x is nan (avoid invalid except in fistp) fld %st(0) fistpl -16(%rsp) fildl -16(%rsp) @@ -58,28 +59,25 @@ exp2l: add %eax,-16(%rsp) f2xm1 fld1 - faddp - fldt -24(%rsp) + faddp # 2^(x-rint(x)) + fldt -24(%rsp) # 2^rint(x) fmulp - fstp %st(1) +2: fstp %st(1) ret -2: fld %st(0) +3: fld %st(0) fstpt -24(%rsp) - mov -15(%rsp),%ah - and $0x7f,%ah - cmp $0x7f,%ah - jne 1f - decb -15(%rsp) - fstp %st(0) - fldt -24(%rsp) -1: fld %st(0) + fld1 + mov -15(%rsp),%ax + and $0x7fff,%ax + cmp $0x7fff,%ax + je 1f # x = +-inf + fld %st(1) frndint - fxch %st(1) - fsub %st(1) + fxch %st(2) + fsub %st(2) # st(0)=x-rint(x), st(1)=1, st(2)=rint(x) f2xm1 - fld1 - faddp - fscale + faddp # 2^(x-rint(x)) +1: fscale fstp %st(1) ret |