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