13a8617a8SJordan K. Hubbard /*
23a8617a8SJordan K. Hubbard * ====================================================
33a8617a8SJordan K. Hubbard * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
43a8617a8SJordan K. Hubbard *
53a8617a8SJordan K. Hubbard * Developed at SunPro, a Sun Microsystems, Inc. business.
63a8617a8SJordan K. Hubbard * Permission to use, copy, modify, and distribute this
73a8617a8SJordan K. Hubbard * software is freely granted, provided that this notice
83a8617a8SJordan K. Hubbard * is preserved.
93a8617a8SJordan K. Hubbard * ====================================================
103a8617a8SJordan K. Hubbard */
113a8617a8SJordan K. Hubbard
123a8617a8SJordan K. Hubbard /* Tanh(x)
133a8617a8SJordan K. Hubbard * Return the Hyperbolic Tangent of x
143a8617a8SJordan K. Hubbard *
153a8617a8SJordan K. Hubbard * Method :
163a8617a8SJordan K. Hubbard * x -x
173a8617a8SJordan K. Hubbard * e - e
183a8617a8SJordan K. Hubbard * 0. tanh(x) is defined to be -----------
193a8617a8SJordan K. Hubbard * x -x
203a8617a8SJordan K. Hubbard * e + e
213a8617a8SJordan K. Hubbard * 1. reduce x to non-negative by tanh(-x) = -tanh(x).
22fe72622eSBruce Evans * 2. 0 <= x < 2**-28 : tanh(x) := x with inexact if x != 0
233a8617a8SJordan K. Hubbard * -t
24fe72622eSBruce Evans * 2**-28 <= x < 1 : tanh(x) := -----; t = expm1(-2x)
253a8617a8SJordan K. Hubbard * t + 2
263a8617a8SJordan K. Hubbard * 2
27fe72622eSBruce Evans * 1 <= x < 22 : tanh(x) := 1 - -----; t = expm1(2x)
283a8617a8SJordan K. Hubbard * t + 2
29fe72622eSBruce Evans * 22 <= x <= INF : tanh(x) := 1.
303a8617a8SJordan K. Hubbard *
313a8617a8SJordan K. Hubbard * Special cases:
323a8617a8SJordan K. Hubbard * tanh(NaN) is NaN;
333a8617a8SJordan K. Hubbard * only tanh(0)=0 is exact for finite argument.
343a8617a8SJordan K. Hubbard */
353a8617a8SJordan K. Hubbard
36a48e1f22SSteve Kargl #include <float.h>
37a48e1f22SSteve Kargl
383a8617a8SJordan K. Hubbard #include "math.h"
393a8617a8SJordan K. Hubbard #include "math_private.h"
403a8617a8SJordan K. Hubbard
41*53c7f228SSteve Kargl static const volatile double tiny = 1.0e-300;
4279c2e69aSSteve Kargl static const double one = 1.0, two = 2.0, huge = 1.0e300;
433a8617a8SJordan K. Hubbard
4459b19ff1SAlfred Perlstein double
tanh(double x)4559b19ff1SAlfred Perlstein tanh(double x)
463a8617a8SJordan K. Hubbard {
473a8617a8SJordan K. Hubbard double t,z;
483a8617a8SJordan K. Hubbard int32_t jx,ix;
493a8617a8SJordan K. Hubbard
503a8617a8SJordan K. Hubbard GET_HIGH_WORD(jx,x);
513a8617a8SJordan K. Hubbard ix = jx&0x7fffffff;
523a8617a8SJordan K. Hubbard
533a8617a8SJordan K. Hubbard /* x is INF or NaN */
543a8617a8SJordan K. Hubbard if(ix>=0x7ff00000) {
553a8617a8SJordan K. Hubbard if (jx>=0) return one/x+one; /* tanh(+-inf)=+-1 */
563a8617a8SJordan K. Hubbard else return one/x-one; /* tanh(NaN) = NaN */
573a8617a8SJordan K. Hubbard }
583a8617a8SJordan K. Hubbard
593a8617a8SJordan K. Hubbard /* |x| < 22 */
603a8617a8SJordan K. Hubbard if (ix < 0x40360000) { /* |x|<22 */
61fe72622eSBruce Evans if (ix<0x3e300000) { /* |x|<2**-28 */
62fe72622eSBruce Evans if(huge+x>one) return x; /* tanh(tiny) = tiny with inexact */
63fe72622eSBruce Evans }
643a8617a8SJordan K. Hubbard if (ix>=0x3ff00000) { /* |x|>=1 */
653a8617a8SJordan K. Hubbard t = expm1(two*fabs(x));
663a8617a8SJordan K. Hubbard z = one - two/(t+two);
673a8617a8SJordan K. Hubbard } else {
683a8617a8SJordan K. Hubbard t = expm1(-two*fabs(x));
693a8617a8SJordan K. Hubbard z= -t/(t+two);
703a8617a8SJordan K. Hubbard }
71fe72622eSBruce Evans /* |x| >= 22, return +-1 */
723a8617a8SJordan K. Hubbard } else {
73fe72622eSBruce Evans z = one - tiny; /* raise inexact flag */
743a8617a8SJordan K. Hubbard }
753a8617a8SJordan K. Hubbard return (jx>=0)? z: -z;
763a8617a8SJordan K. Hubbard }
77a48e1f22SSteve Kargl
78a48e1f22SSteve Kargl #if (LDBL_MANT_DIG == 53)
79a48e1f22SSteve Kargl __weak_reference(tanh, tanhl);
80a48e1f22SSteve Kargl #endif
81