fluencelabs/musl
1
0
Fork 0
mirror of https://github.com/fluencelabs/musl synced 2025-06-28 22:22:01 +00:00
Code Issues Projects Releases Wiki Activity

math: fix x86 expl.s to raise underflow and clean up special case handling

Browse Source
This commit is contained in:
Szabolcs Nagy
2013-08-15 10:54:56 +00:00
parent c599f4f4d7
commit 1b3973fb43
2 changed files with 31 additions and 45 deletions
Download Patch File Download Diff File Expand all files Collapse all files

39
src/math/i386/expl.s
View File

@ -8,34 +8,27 @@
expl: expl:
fldt 4(%esp) fldt 4(%esp)
# special cases: 2*x is +-inf, nan or |x| < 0x1p-32 # interesting case: 0x1p-32 <= |x| < 16384
# check (exponent|0x8000)+2 < 0xbfff+2-32 # check if (exponent|0x8000) is in [0xbfff-32, 0xbfff+13]
movw 12(%esp), %ax mov 12(%esp), %ax
movw %ax, %dx or $0x8000, %ax
orw $0x8000, %dx sub $0xbfdf, %ax
addw $2, %dx cmp $45, %ax
cmpw $0xbfff-30, %dx jbe 2f
jnb 3f test %ax, %ax
cmpw $1, %dx
jbe 1f
# if |x|<0x1p-32 return 1+x
fld1 fld1
jmp 2f js 1f
1: testw %ax, %ax # if |x|>=0x1p14 or nan return 2^trunc(x)
jns 1f fscale
# if 2*x == -inf,-nan return -0/x fstp %st(1)
fldz
fchs
fdivp
ret ret
# if 2*x == inf,nan return 2*x # if |x|<0x1p-32 return 1+x
1: fld %st(0) 1: faddp
2: faddp
ret ret
# should be 0x1.71547652b82fe178p0 == 0x3fff b8aa3b29 5c17f0bc # should be 0x1.71547652b82fe178p0L == 0x3fff b8aa3b29 5c17f0bc
# it will be wrong on non-nearest rounding mode # it will be wrong on non-nearest rounding mode
3: fldl2e 2: fldl2e
subl $44, %esp subl $44, %esp
# hi = log2e_hi*x # hi = log2e_hi*x
# 2^hi = exp2l(hi) # 2^hi = exp2l(hi)

37
src/math/x86_64/expl.s
View File

@ -8,32 +8,25 @@
expl: expl:
fldt 8(%rsp) fldt 8(%rsp)
# special cases: 2*x is +-inf, nan or |x| < 0x1p-32 # interesting case: 0x1p-32 <= |x| < 16384
# check (exponent|0x8000)+2 < 0xbfff+2-32 # check if (exponent|0x8000) is in [0xbfff-32, 0xbfff+13]
movw 16(%rsp), %ax mov 16(%rsp), %ax
movw %ax, %dx or $0x8000, %ax
orw $0x8000, %dx sub $0xbfdf, %ax
addw $2, %dx cmp $45, %ax
cmpw $0xbfff-30, %dx jbe 2f
jnb 3f test %ax, %ax
cmpw $1, %dx
jbe 1f
# if |x|<0x1p-32 return 1+x
fld1 fld1
jmp 2f js 1f
1: testw %ax, %ax # if |x|>=0x1p14 or nan return 2^trunc(x)
jns 1f fscale
# if 2*x == -inf,-nan return -0/x fstp %st(1)
fldz
fchs
fdivp
ret ret
# if 2*x == inf,nan return 2*x # if |x|<0x1p-32 return 1+x
1: fld %st(0) 1: faddp
2: faddp
ret ret
# should be 0x1.71547652b82fe178p0 == 0x3fff b8aa3b29 5c17f0bc # should be 0x1.71547652b82fe178p0L == 0x3fff b8aa3b29 5c17f0bc
# it will be wrong on non-nearest rounding mode # it will be wrong on non-nearest rounding mode
3: fldl2e 3: fldl2e
subq $48, %rsp subq $48, %rsp