ATANH

[Damian McGuckin <damianm@esi.com.au>]

The comments for this routine say:

  * atanh(x) = log((1+x)/(1-x))/2 = log1p(2x/(1-x))/2 ~= x + x^3/3 + o(x^5)

There is a point where atanh(x) can be approximated by just x to machine 
precision. This is where the exponent of x is less than some number, or 
where x itself is less than some number.

In MUSL, for the 80-bit version, directly from the code, this is

 	2^(-LDBL_MANT_DIG/2) = 0x1.0p-32;

Interestingly, this same number, is  used double atanh and float atanhf.

Note that in fdlibm,  0x1.0p-28 is used for all types

In fact, isn't

 		x * (1 + x^2/3) == x (after roundoff)

if
 	x < 2^(p/2),

i.e. 0x1.0p-26 for double, 0x1.0p-12 for float.

Note that for 80-bit wide, it is 0x1.0p-32 agreeing with MUSL!!

Regards - Damian

Pacific Engineering Systems International, 277-279 Broadway, Glebe NSW 2037
Ph:+61-2-8571-0847 .. Fx:+61-2-9692-9623 | unsolicited email not wanted here
Views & opinions here are mine and not those of any past or present employer