Category : C Source Code
Archive   : CEPHES22.ZIP
Filename : ATANHL.C

/* atanhl.c
*
* Inverse hyperbolic tangent, long double precision
*
*
*
* SYNOPSIS:
*
* long double x, y, atanhl();
*
* y = atanhl( x );
*
*
*
* DESCRIPTION:
*
* Returns inverse hyperbolic tangent of argument in the range
* MINLOGL to MAXLOGL.
*
* If |x| < 0.5, the rational form x + x**3 P(x)/Q(x) is
* employed. Otherwise,
* atanh(x) = 0.5 * log( (1+x)/(1-x) ).
*
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* IEEE -1,1 30000 1.1e-19 3.3e-20
*
*/



/*
Cephes Math Library Release 2.2: January, 1991
Copyright (C) 1987, 1991 by Stephen L. Moshier
Direct inquiries to 30 Frost Street, Cambridge, MA 02140
*/

#include "mconf.h"

#ifdef UNK
static long double P[] = {
2.9647757819596835680719E-3L,
-8.0026596513099094380633E-1L,
7.7920941408493040219831E0L,
-2.4330686602187898836837E1L,
3.0204265014595622991082E1L,
-1.2961142942114056581210E1L,
};
static long double Q[] = {
/* 1.0000000000000000000000E0L,*/
-1.3729634163247557081869E1L,
6.2320841104088512332185E1L,
-1.2469344457045341444078E2L,
1.1394285233959210574352E2L,
-3.8883428826342169425890E1L,
};
#endif

#ifdef IBMPC
static short P[] = {
0x3aa2,0x036b,0xaf06,0xc24c,0x3ff6,
0x528e,0x56e8,0x3af4,0xccde,0xbffe,
0x9d89,0xc9a1,0xd5cf,0xf958,0x4001,
0xa653,0x6cfa,0x3f04,0xc2a5,0xc003,
0xc651,0x2b3d,0x55b2,0xf1a2,0x4003,
0xd76d,0xf293,0xd76b,0xcf60,0xc002,
};
static short Q[] = {
/*0x0000,0x0000,0x0000,0x8000,0x3fff,*/
0xd1b9,0x5314,0x94df,0xdbac,0xc002,
0x3caa,0x0517,0x8a92,0xf948,0x4004,
0x535e,0xaf5f,0x0b2a,0xf963,0xc005,
0xa6f9,0xb702,0xbd8a,0xe3e2,0x4005,
0xe136,0xf5ee,0xa190,0x9b88,0xc004,
};
#endif

#ifdef MIEEE
static long P[] = {
0x3ff60000,0xc24caf06,0x036b3aa2,
0xbffe0000,0xccde3af4,0x56e8528e,
0x40010000,0xf958d5cf,0xc9a19d89,
0xc0030000,0xc2a53f04,0x6cfaa653,
0x40030000,0xf1a255b2,0x2b3dc651,
0xc0020000,0xcf60d76b,0xf293d76d,
};
static long Q[] = {
/*0x3fff0000,0x80000000,0x00000000,*/
0xc0020000,0xdbac94df,0x5314d1b9,
0x40040000,0xf9488a92,0x05173caa,
0xc0050000,0xf9630b2a,0xaf5f535e,
0x40050000,0xe3e2bd8a,0xb702a6f9,
0xc0040000,0x9b88a190,0xf5eee136,
};
#endif

extern double MAXNUML;

long double atanhl(x)
long double x;
{
long double s, z;
int i, n;
long double fabsl(), logl(), polevll(), p1evll();

z = fabsl(x);
if( z >= 1.0L )
{
if( x == 1.0L )
return( MAXNUML );
if( x == -1.0L )
return( -MAXNUML );
mtherr( "atanhl", DOMAIN );
return( MAXNUML );
}

if( z < 1.0e-8L )
return(x);

if( z < 0.5L )
{
z = x * x;
s = x + x * z * (polevll(z, P, 5) / p1evll(z, Q, 5));
return(s);
}

return( 0.5L * logl((1.0L+x)/(1.0L-x)) );
}