diff options
| author | Szabolcs Nagy <nsz@port70.net> | 2013-09-03 15:02:10 +0000 |
|---|---|---|
| committer | Szabolcs Nagy <nsz@port70.net> | 2013-09-05 11:30:07 +0000 |
| commit | bcd797a5ba4631c031919dad832d670e564212e9 (patch) | |
| tree | 4ad873585d3e87dfca85b16baa3b6613a7e9f479 /src/math/atan2l.c | |
| parent | c2a0dfea629617a39af2f59bd400e1a3595d0783 (diff) | |
| download | musl-bcd797a5ba4631c031919dad832d670e564212e9.tar.gz musl-bcd797a5ba4631c031919dad832d670e564212e9.tar.xz musl-bcd797a5ba4631c031919dad832d670e564212e9.zip | |
math: long double inverse trigonometric cleanup (acosl, asinl, atanl, atan2l)
* added ld128 support from freebsd fdlibm (untested) * using new ldshape union instead of IEEEl2bits * inexact status flag is not supported
Diffstat (limited to 'src/math/atan2l.c')
| -rw-r--r-- | src/math/atan2l.c | 50 |
1 files changed, 22 insertions, 28 deletions
diff --git a/src/math/atan2l.c b/src/math/atan2l.c index e0167d09..f0937a97 100644 --- a/src/math/atan2l.c +++ b/src/math/atan2l.c @@ -27,56 +27,50 @@ long double atan2l(long double y, long double x) long double atan2l(long double y, long double x) { - union IEEEl2bits ux, uy; + union ldshape ux, uy; long double z; - int m; - uint16_t exptx, expsignx, expty, expsigny; + int m, ex, ey; if (isnan(x) || isnan(y)) return x+y; if (x == 1) return atanl(y); - uy.e = y; - expsigny = uy.xbits.expsign; - expty = expsigny & 0x7fff; - ux.e = x; - expsignx = ux.xbits.expsign; - exptx = expsignx & 0x7fff; - m = ((expsigny>>15)&1) | ((expsignx>>14)&2); /* 2*sign(x)+sign(y) */ - + ux.f = x; + uy.f = y; + ex = ux.i.se & 0x7fff; + ey = uy.i.se & 0x7fff; + m = 2*(ux.i.se>>15) | uy.i.se>>15; if (y == 0) { switch(m) { case 0: case 1: return y; /* atan(+-0,+anything)=+-0 */ - case 2: return 2*pio2_hi+0x1p-120f; /* atan(+0,-anything) = pi */ - case 3: return -2*pio2_hi-0x1p-120f; /* atan(-0,-anything) =-pi */ + case 2: return 2*pio2_hi; /* atan(+0,-anything) = pi */ + case 3: return -2*pio2_hi; /* atan(-0,-anything) =-pi */ } } if (x == 0) - return m&1 ? -pio2_hi-0x1p-120f : pio2_hi+0x1p-120f; - /* when x is INF */ - if (exptx == 0x7fff) { - if (expty == 0x7fff) { + return m&1 ? -pio2_hi : pio2_hi; + if (ex == 0x7fff) { + if (ey == 0x7fff) { switch(m) { - case 0: return pio2_hi*0.5+0x1p-120f; /* atan(+INF,+INF) */ - case 1: return -pio2_hi*0.5-0x1p-120f; /* atan(-INF,+INF) */ - case 2: return 1.5*pio2_hi+0x1p-120f; /* atan(+INF,-INF) */ - case 3: return -1.5*pio2_hi-0x1p-120f; /* atan(-INF,-INF) */ + case 0: return pio2_hi/2; /* atan(+INF,+INF) */ + case 1: return -pio2_hi/2; /* atan(-INF,+INF) */ + case 2: return 1.5*pio2_hi; /* atan(+INF,-INF) */ + case 3: return -1.5*pio2_hi; /* atan(-INF,-INF) */ } } else { switch(m) { case 0: return 0.0; /* atan(+...,+INF) */ case 1: return -0.0; /* atan(-...,+INF) */ - case 2: return 2*pio2_hi+0x1p-120f; /* atan(+...,-INF) */ - case 3: return -2*pio2_hi-0x1p-120f; /* atan(-...,-INF) */ + case 2: return 2*pio2_hi; /* atan(+...,-INF) */ + case 3: return -2*pio2_hi; /* atan(-...,-INF) */ } } } - /* when y is INF */ - if (exptx+120 < expty || expty == 0x7fff) - return m&1 ? -pio2_hi-0x1p-120f : pio2_hi+0x1p-120f; - - if ((m&2) && expty+120 < exptx) /* |y/x| tiny, x<0 */ + if (ex+120 < ey || ey == 0x7fff) + return m&1 ? -pio2_hi : pio2_hi; + /* z = atan(|y/x|) without spurious underflow */ + if ((m&2) && ey+120 < ex) /* |y/x| < 0x1p-120, x<0 */ z = 0.0; else z = atanl(fabsl(y/x)); |