xref: /csrg-svn/old/libm/libm/tanh.c (revision 21035)
1*21035Smiriam /*
2*21035Smiriam  * Copyright (c) 1985 Regents of the University of California.
3*21035Smiriam  *
4*21035Smiriam  * Use and reproduction of this software are granted  in  accordance  with
5*21035Smiriam  * the terms and conditions specified in  the  Berkeley  Software  License
6*21035Smiriam  * Agreement (in particular, this entails acknowledgement of the programs'
7*21035Smiriam  * source, and inclusion of this notice) with the additional understanding
8*21035Smiriam  * that  all  recipients  should regard themselves as participants  in  an
9*21035Smiriam  * ongoing  research  project and hence should  feel  obligated  to report
10*21035Smiriam  * their  experiences (good or bad) with these elementary function  codes,
11*21035Smiriam  * using "sendbug 4bsd-bugs@BERKELEY", to the authors.
12*21035Smiriam  */
139942Ssam 
14*21035Smiriam #ifndef lint
15*21035Smiriam static char sccsid[] = "@(#)tanh.c	4.2 (Berkeley) 05/23/85";
16*21035Smiriam #endif not lint
179942Ssam 
18*21035Smiriam /* TANH(X)
19*21035Smiriam  * RETURN THE HYPERBOLIC TANGENT OF X
20*21035Smiriam  * DOUBLE PRECISION (VAX D FORMAT 56 BITS, IEEE DOUBLE 53 BITS)
21*21035Smiriam  * CODED IN C BY K.C. NG, 1/8/85;
22*21035Smiriam  * REVISED BY K.C. NG on 2/8/85, 2/11/85, 3/7/85, 3/24/85.
23*21035Smiriam  *
24*21035Smiriam  * Required system supported functions :
25*21035Smiriam  *	copysign(x,y)
26*21035Smiriam  *	finite(x)
27*21035Smiriam  *
28*21035Smiriam  * Required kernel function:
29*21035Smiriam  *	E(x)	...exp(x)-1
30*21035Smiriam  *
31*21035Smiriam  * Method :
32*21035Smiriam  *	1. reduce x to non-negative by tanh(-x) = - tanh(x).
33*21035Smiriam  *	2. For appropriate values of small,
34*21035Smiriam  *					          -E(-2x)
35*21035Smiriam  *	    0     <  x <=     1   :  tanh(x) := ------------
36*21035Smiriam  *					         E(-2x) + 2
37*21035Smiriam  *							 2
38*21035Smiriam  *	    1     <= x <= 22.0    :  tanh(x) := 1 -  ------------
39*21035Smiriam  *						      E(2x) + 2
40*21035Smiriam  *	    22.0  <  x <= INF     :  tanh(x) := 1.
41*21035Smiriam  *
42*21035Smiriam  *	Note: 22 are chosen so that fl(1.0+2/(E(2*22)+2)) == 1.
43*21035Smiriam  *
44*21035Smiriam  * Special cases:
45*21035Smiriam  *	tanh(NAN) is NAN;
46*21035Smiriam  *	only tanh(0)=0 is exact for finite argument.
47*21035Smiriam  *
48*21035Smiriam  * Accuracy:
49*21035Smiriam  *	tanh(x) returns the exact hyperbolic tangent of x nealy rounded.
50*21035Smiriam  *	In a test run with 1,024,000 random arguments on a VAX, the maximum
51*21035Smiriam  *	observed error was 2.22 ulps (units in the last place).
52*21035Smiriam  */
539942Ssam 
54*21035Smiriam double tanh(x)
55*21035Smiriam double x;
569942Ssam {
57*21035Smiriam 	static double one=1.0, two=2.0;
58*21035Smiriam 	double E(), t,  copysign(), sign;
59*21035Smiriam 	int finite();
609942Ssam 
61*21035Smiriam 	if(x!=x) return(x);
629942Ssam 
63*21035Smiriam 	sign=copysign(one,x);
64*21035Smiriam 	x=copysign(x,one);
65*21035Smiriam 	if(x < 22.0)
66*21035Smiriam 	    if( x > one )
67*21035Smiriam 		return(copysign(one-two/(E(x+x)+two),sign));
68*21035Smiriam 	    else
69*21035Smiriam 		{t= -E(-(x+x)); return(copysign(t/(two-t),sign));}
709942Ssam 
71*21035Smiriam 	else if(finite(x))
72*21035Smiriam 	    return (sign+1.0E-37); /* raise the inexact flag */
73*21035Smiriam 
74*21035Smiriam 	else
75*21035Smiriam 	    return(sign);	/* x is +- INF */
769942Ssam }
77*21035Smiriam 
78