X-Git-Url: http://nsz.repo.hu/git/?a=blobdiff_plain;f=src%2Fmath%2Fasinl.c;h=347c53568cb11d787f91caaf21b2673e57b29da8;hb=d055e6a45a17673b8dd3ec16e786bb2fe1700dd5;hp=8799341d568fd81a8d94d439b3c86297fe644566;hpb=5d5ab51862cbd010bdf52dc3b04b0967450bcd1a;p=musl diff --git a/src/math/asinl.c b/src/math/asinl.c index 8799341d..347c5356 100644 --- a/src/math/asinl.c +++ b/src/math/asinl.c @@ -23,27 +23,29 @@ long double asinl(long double x) } #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384 #include "__invtrigl.h" -/* 0.95 */ -#define THRESH ((0xe666666666666666ULL>>(64-(LDBL_MANH_SIZE-1)))|LDBL_NBIT) +#if LDBL_MANT_DIG == 64 +#define CLOSETO1(u) (u.i.m>>56 >= 0xf7) +#define CLEARBOTTOM(u) (u.i.m &= -1ULL << 32) +#elif LDBL_MANT_DIG == 113 +#define CLOSETO1(u) (u.i.top >= 0xee00) +#define CLEARBOTTOM(u) (u.i.lo = 0) +#endif long double asinl(long double x) { - union IEEEl2bits u; - long double z,r,s; - uint16_t expsign, expt; + union ldshape u = {x}; + long double z, r, s; + uint16_t e = u.i.se & 0x7fff; + int sign = u.i.se >> 15; - u.e = x; - expsign = u.xbits.expsign; - expt = expsign & 0x7fff; - if (expt >= 0x3fff) { /* |x| >= 1 or nan */ - if (expt == 0x3fff && - ((u.bits.manh&~LDBL_NBIT)|u.bits.manl) == 0) - /* asin(+-1)=+-pi/2 with inexact */ + if (e >= 0x3fff) { /* |x| >= 1 or nan */ + /* asin(+-1)=+-pi/2 with inexact */ + if (x == 1 || x == -1) return x*pio2_hi + 0x1p-120f; return 0/(x-x); } - if (expt < 0x3fff - 1) { /* |x| < 0.5 */ - if (expt < 0x3fff - 32) { /* |x|<0x1p-32, asinl(x)=x */ + if (e < 0x3fff - 1) { /* |x| < 0.5 */ + if (e < 0x3fff - (LDBL_MANT_DIG+1)/2) { /* return x with inexact if x!=0 */ FORCE_EVAL(x + 0x1p120f); return x; @@ -54,18 +56,16 @@ long double asinl(long double x) z = (1.0 - fabsl(x))*0.5; s = sqrtl(z); r = __invtrigl_R(z); - if (u.bits.manh >= THRESH) { /* if |x| is close to 1 */ + if (CLOSETO1(u)) { x = pio2_hi - (2*(s+s*r)-pio2_lo); } else { long double f, c; - u.e = s; - u.bits.manl = 0; - f = u.e; - c = (z-f*f)/(s+f); + u.f = s; + CLEARBOTTOM(u); + f = u.f; + c = (z - f*f)/(s + f); x = 0.5*pio2_hi-(2*s*r - (pio2_lo-2*c) - (0.5*pio2_hi-2*f)); } - if (expsign>>15) - return -x; - return x; + return sign ? -x : x; } #endif