1#include "libm.h" 2 3#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 4long double atanhl(long double x) { 5 return atanh(x); 6} 7#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384 8/* atanh(x) = log((1+x)/(1-x))/2 = log1p(2x/(1-x))/2 ~= x + x^3/3 + o(x^5) */ 9long double atanhl(long double x) { 10 union ldshape u = {x}; 11 unsigned e = u.i.se & 0x7fff; 12 unsigned s = u.i.se >> 15; 13 14 /* |x| */ 15 u.i.se = e; 16 x = u.f; 17 18 if (e < 0x3ff - 1) { 19 if (e < 0x3ff - LDBL_MANT_DIG / 2) { 20 /* handle underflow */ 21 if (e == 0) 22 FORCE_EVAL((float)x); 23 } else { 24 /* |x| < 0.5, up to 1.7ulp error */ 25 x = 0.5 * log1pl(2 * x + 2 * x * x / (1 - x)); 26 } 27 } else { 28 /* avoid overflow */ 29 x = 0.5 * log1pl(2 * (x / (1 - x))); 30 } 31 return s ? -x : x; 32} 33#endif 34