xref: /illumos-gate/usr/src/lib/libm/common/Q/tanhl.c (revision 5422785d352a2bb398daceab3d1898a8aa64d006)
1 /*
2  * CDDL HEADER START
3  *
4  * The contents of this file are subject to the terms of the
5  * Common Development and Distribution License (the "License").
6  * You may not use this file except in compliance with the License.
7  *
8  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9  * or http://www.opensolaris.org/os/licensing.
10  * See the License for the specific language governing permissions
11  * and limitations under the License.
12  *
13  * When distributing Covered Code, include this CDDL HEADER in each
14  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15  * If applicable, add the following below this CDDL HEADER, with the
16  * fields enclosed by brackets "[]" replaced with your own identifying
17  * information: Portions Copyright [yyyy] [name of copyright owner]
18  *
19  * CDDL HEADER END
20  */
21 
22 /*
23  * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
24  */
25 /*
26  * Copyright 2006 Sun Microsystems, Inc.  All rights reserved.
27  * Use is subject to license terms.
28  */
29 
30 #pragma weak __tanhl = tanhl
31 
32 /*
33  * tanhl(x) returns the hyperbolic tangent of x
34  *
35  * Method :
36  *	1. reduce x to non-negative:  tanhl(-x) = - tanhl(x).
37  *	2.
38  *	  0      <  x <=  small    :  tanhl(x) := x
39  *					          -expm1l(-2x)
40  *	  small  <  x <=  1        :  tanhl(x) := --------------
41  *					         expm1l(-2x) + 2
42  *							  2
43  *	  1      <= x <= threshold :  tanhl(x) := 1 -  ---------------
44  *						      expm1l(2x) + 2
45  *     threshold <  x <= INF       :  tanhl(x) := 1.
46  *
47  * where
48  *	single : 	small = 1.e-5		threshold = 11.0
49  *	double : 	small = 1.e-10		threshold = 22.0
50  *	quad   : 	small = 1.e-20		threshold = 45.0
51  *
52  * Note: threshold was chosen so that
53  *		fl(1.0+2/(expm1(2*threshold)+2)) == 1.
54  *
55  * Special cases:
56  *	tanhl(NaN) is NaN;
57  *	only tanhl(0.0)=0.0 is exact for finite argument.
58  */
59 
60 #include "libm.h"
61 #include "longdouble.h"
62 
63 static const long double small = 1.0e-20L, one = 1.0, two = 2.0,
64 #ifndef lint
65 	big = 1.0e+20L,
66 #endif
67 	threshold = 45.0L;
68 
69 long double
70 tanhl(long double x) {
71 	long double t, y, z;
72 	int signx;
73 	volatile long double dummy;
74 
75 	if (isnanl(x))
76 		return (x + x);		/* x is NaN */
77 	signx = signbitl(x);
78 	t = fabsl(x);
79 	z = one;
80 	if (t <= threshold) {
81 		if (t > one)
82 			z = one - two / (expm1l(t + t) + two);
83 		else if (t > small) {
84 			y = expm1l(-t - t);
85 			z = -y / (y + two);
86 		} else {
87 #ifndef lint
88 			dummy = t + big;
89 							/* inexact if t != 0 */
90 #endif
91 			return (x);
92 		}
93 	} else if (!finitel(t))
94 		return (copysignl(one, x));
95 	else
96 		return (signx ? -z + small * small : z - small * small);
97 	return (signx ? -z : z);
98 }
99