1*627f7eb2Smrg /* s_asinhl.c -- long double version of s_asinh.c.
2*627f7eb2Smrg * Conversion to long double by Ulrich Drepper,
3*627f7eb2Smrg * Cygnus Support, drepper@cygnus.com.
4*627f7eb2Smrg */
5*627f7eb2Smrg
6*627f7eb2Smrg /*
7*627f7eb2Smrg * ====================================================
8*627f7eb2Smrg * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
9*627f7eb2Smrg *
10*627f7eb2Smrg * Developed at SunPro, a Sun Microsystems, Inc. business.
11*627f7eb2Smrg * Permission to use, copy, modify, and distribute this
12*627f7eb2Smrg * software is freely granted, provided that this notice
13*627f7eb2Smrg * is preserved.
14*627f7eb2Smrg * ====================================================
15*627f7eb2Smrg */
16*627f7eb2Smrg
17*627f7eb2Smrg #if defined(LIBM_SCCS) && !defined(lint)
18*627f7eb2Smrg static char rcsid[] = "NetBSD: ";
19*627f7eb2Smrg #endif
20*627f7eb2Smrg
21*627f7eb2Smrg /* asinhq(x)
22*627f7eb2Smrg * Method :
23*627f7eb2Smrg * Based on
24*627f7eb2Smrg * asinhq(x) = signl(x) * logq [ |x| + sqrtq(x*x+1) ]
25*627f7eb2Smrg * we have
26*627f7eb2Smrg * asinhq(x) := x if 1+x*x=1,
27*627f7eb2Smrg * := signl(x)*(logq(x)+ln2)) for large |x|, else
28*627f7eb2Smrg * := signl(x)*logq(2|x|+1/(|x|+sqrtq(x*x+1))) if|x|>2, else
29*627f7eb2Smrg * := signl(x)*log1pq(|x| + x^2/(1 + sqrtq(1+x^2)))
30*627f7eb2Smrg */
31*627f7eb2Smrg
32*627f7eb2Smrg #include "quadmath-imp.h"
33*627f7eb2Smrg
34*627f7eb2Smrg static const __float128
35*627f7eb2Smrg one = 1,
36*627f7eb2Smrg ln2 = 6.931471805599453094172321214581765681e-1Q,
37*627f7eb2Smrg huge = 1.0e+4900Q;
38*627f7eb2Smrg
39*627f7eb2Smrg __float128
asinhq(__float128 x)40*627f7eb2Smrg asinhq (__float128 x)
41*627f7eb2Smrg {
42*627f7eb2Smrg __float128 t, w;
43*627f7eb2Smrg int32_t ix, sign;
44*627f7eb2Smrg ieee854_float128 u;
45*627f7eb2Smrg
46*627f7eb2Smrg u.value = x;
47*627f7eb2Smrg sign = u.words32.w0;
48*627f7eb2Smrg ix = sign & 0x7fffffff;
49*627f7eb2Smrg if (ix == 0x7fff0000)
50*627f7eb2Smrg return x + x; /* x is inf or NaN */
51*627f7eb2Smrg if (ix < 0x3fc70000)
52*627f7eb2Smrg { /* |x| < 2^ -56 */
53*627f7eb2Smrg math_check_force_underflow (x);
54*627f7eb2Smrg if (huge + x > one)
55*627f7eb2Smrg return x; /* return x inexact except 0 */
56*627f7eb2Smrg }
57*627f7eb2Smrg u.words32.w0 = ix;
58*627f7eb2Smrg if (ix > 0x40350000)
59*627f7eb2Smrg { /* |x| > 2 ^ 54 */
60*627f7eb2Smrg w = logq (u.value) + ln2;
61*627f7eb2Smrg }
62*627f7eb2Smrg else if (ix >0x40000000)
63*627f7eb2Smrg { /* 2^ 54 > |x| > 2.0 */
64*627f7eb2Smrg t = u.value;
65*627f7eb2Smrg w = logq (2.0 * t + one / (sqrtq (x * x + one) + t));
66*627f7eb2Smrg }
67*627f7eb2Smrg else
68*627f7eb2Smrg { /* 2.0 > |x| > 2 ^ -56 */
69*627f7eb2Smrg t = x * x;
70*627f7eb2Smrg w = log1pq (u.value + t / (one + sqrtq (one + t)));
71*627f7eb2Smrg }
72*627f7eb2Smrg if (sign & 0x80000000)
73*627f7eb2Smrg return -w;
74*627f7eb2Smrg else
75*627f7eb2Smrg return w;
76*627f7eb2Smrg }
77