xref: /illumos-gate/usr/src/lib/libm/common/C/tanh.c (revision e4fb8a5fd5bb5a5d316430c6ac8f6d33704fcf5d)
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 __tanh = tanh
31 
32 /* INDENT OFF */
33 /*
34  * TANH(X)
35  * RETURN THE HYPERBOLIC TANGENT OF X
36  * code based on 4.3bsd
37  * Modified by K.C. Ng for sun 4.0, Jan 31, 1987
38  *
39  * Method :
40  *	1. reduce x to non-negative by tanh(-x) = - tanh(x).
41  *	2.
42  *	    0      <  x <=  1.e-10 :  tanh(x) := x
43  *					          -expm1(-2x)
44  *	    1.e-10 <  x <=  1      :  tanh(x) := --------------
45  *					         expm1(-2x) + 2
46  *							  2
47  *	    1      <= x <=  22.0   :  tanh(x) := 1 -  ---------------
48  *						      expm1(2x) + 2
49  *	    22.0   <  x <= INF     :  tanh(x) := 1.
50  *
51  *	Note: 22 was chosen so that fl(1.0+2/(expm1(2*22)+2)) == 1.
52  *
53  * Special cases:
54  *	tanh(NaN) is NaN;
55  *	only tanh(0)=0 is exact for finite argument.
56  */
57 
58 #include "libm.h"
59 #include "libm_protos.h"
60 #include <math.h>
61 
62 static const double
63 	one = 1.0,
64 	two = 2.0,
65 	small = 1.0e-10,
66 	big = 1.0e10;
67 /* INDENT ON */
68 
69 double
70 tanh(double x)
71 {
72 	double t, y, z;
73 	int signx;
74 	volatile double dummy __unused;
75 
76 	if (isnan(x))
77 		return (x * x);	/* + -> * for Cheetah */
78 	signx = signbit(x);
79 	t = fabs(x);
80 	z = one;
81 	if (t <= 22.0) {
82 		if (t > one)
83 			z = one - two / (expm1(t + t) + two);
84 		else if (t > small) {
85 			y = expm1(-t - t);
86 			z = -y / (y + two);
87 		} else {
88 			/* raise the INEXACT flag for non-zero t */
89 			dummy = t + big;
90 #ifdef lint
91 			dummy = dummy;
92 #endif
93 			return (x);
94 		}
95 	} else if (!finite(t))
96 		return (copysign(1.0, x));
97 	else
98 		return ((signx != 0) ? -z + small * small : z - small * small);
99 
100 	return ((signx != 0) ? -z : z);
101 }
102