Back to index

glibc  2.9
e_atanhl.c
Go to the documentation of this file.
00001 /* s_atanhl.c -- long double version of s_atan.c.
00002  * Conversion to long double by Ulrich Drepper,
00003  * Cygnus Support, drepper@cygnus.com.
00004  */
00005 
00006 /*
00007  * ====================================================
00008  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
00009  *
00010  * Developed at SunPro, a Sun Microsystems, Inc. business.
00011  * Permission to use, copy, modify, and distribute this
00012  * software is freely granted, provided that this notice
00013  * is preserved.
00014  * ====================================================
00015  */
00016 
00017 #if defined(LIBM_SCCS) && !defined(lint)
00018 static char rcsid[] = "$NetBSD: $";
00019 #endif
00020 
00021 /* __ieee754_atanhl(x)
00022  * Method :
00023  *    1.Reduced x to positive by atanh(-x) = -atanh(x)
00024  *    2.For x>=0.5
00025  *                   1              2x                          x
00026  *     atanhl(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
00027  *                   2             1 - x                      1 - x
00028  *
00029  *     For x<0.5
00030  *     atanhl(x) = 0.5*log1pl(2x+2x*x/(1-x))
00031  *
00032  * Special cases:
00033  *     atanhl(x) is NaN if |x| > 1 with signal;
00034  *     atanhl(NaN) is that NaN with no signal;
00035  *     atanhl(+-1) is +-INF with signal.
00036  *
00037  */
00038 
00039 #include "math.h"
00040 #include "math_private.h"
00041 
00042 #ifdef __STDC__
00043 static const long double one = 1.0L, huge = 1e4900L;
00044 #else
00045 static long double one = 1.0L, huge = 1e4900L;
00046 #endif
00047 
00048 #ifdef __STDC__
00049 static const long double zero = 0.0L;
00050 #else
00051 static double long zero = 0.0L;
00052 #endif
00053 
00054 #ifdef __STDC__
00055        long double __ieee754_atanhl(long double x)
00056 #else
00057        long double __ieee754_atanhl(x)
00058        long double x;
00059 #endif
00060 {
00061        long double t;
00062        u_int32_t jx, ix;
00063        ieee854_long_double_shape_type u;
00064 
00065        u.value = x;
00066        jx = u.parts32.w0;
00067        ix = jx & 0x7fffffff;
00068        u.parts32.w0 = ix;
00069        if (ix >= 0x3fff0000) /* |x| >= 1.0 or infinity or NaN */
00070          {
00071            if (u.value == one)
00072              return x/zero;
00073            else
00074              return (x-x)/(x-x);
00075          }
00076        if(ix<0x3fc60000 && (huge+x)>zero) return x;     /* x < 2^-57 */
00077 
00078        if(ix<0x3ffe0000) {         /* x < 0.5 */
00079            t = u.value+u.value;
00080            t = 0.5*__log1pl(t+t*u.value/(one-u.value));
00081        } else
00082            t = 0.5*__log1pl((u.value+u.value)/(one-u.value));
00083        if(jx & 0x80000000) return -t; else return t;
00084 }