source:
EcnlProtoTool/trunk/musl-1.1.18/src/math/atanh.c@
444
Last change on this file since 444 was 444, checked in by , 4 years ago | |
---|---|
|
|
File size: 577 bytes |
Rev | Line | |
---|---|---|
[444] | 1 | #include "libm.h" |
2 | ||
3 | /* atanh(x) = log((1+x)/(1-x))/2 = log1p(2x/(1-x))/2 ~= x + x^3/3 + o(x^5) */ | |
4 | double atanh(double x) | |
5 | { | |
6 | union {double f; uint64_t i;} u = {.f = x}; | |
7 | unsigned e = u.i >> 52 & 0x7ff; | |
8 | unsigned s = u.i >> 63; | |
9 | double_t y; | |
10 | ||
11 | /* |x| */ | |
12 | u.i &= (uint64_t)-1/2; | |
13 | y = u.f; | |
14 | ||
15 | if (e < 0x3ff - 1) { | |
16 | if (e < 0x3ff - 32) { | |
17 | /* handle underflow */ | |
18 | if (e == 0) | |
19 | FORCE_EVAL((float)y); | |
20 | } else { | |
21 | /* |x| < 0.5, up to 1.7ulp error */ | |
22 | y = 0.5*log1p(2*y + 2*y*y/(1-y)); | |
23 | } | |
24 | } else { | |
25 | /* avoid overflow */ | |
26 | y = 0.5*log1p(2*(y/(1-y))); | |
27 | } | |
28 | return s ? -y : y; | |
29 | } |
Note:
See TracBrowser
for help on using the repository browser.