1*24590Szliu /* 2*24590Szliu * Copyright (c) 1985 Regents of the University of California. 3*24590Szliu * 4*24590Szliu * Use and reproduction of this software are granted in accordance with 5*24590Szliu * the terms and conditions specified in the Berkeley Software License 6*24590Szliu * Agreement (in particular, this entails acknowledgement of the programs' 7*24590Szliu * source, and inclusion of this notice) with the additional understanding 8*24590Szliu * that all recipients should regard themselves as participants in an 9*24590Szliu * ongoing research project and hence should feel obligated to report 10*24590Szliu * their experiences (good or bad) with these elementary function codes, 11*24590Szliu * using "sendbug 4bsd-bugs@BERKELEY", to the authors. 12*24590Szliu */ 13*24590Szliu 14*24590Szliu #ifndef lint 15*24590Szliu static char sccsid[] = "@(#)atanh.c 1.1 (ELEFUNT) 09/06/85"; 16*24590Szliu #endif not lint 17*24590Szliu 18*24590Szliu /* ATANH(X) 19*24590Szliu * RETURN THE HYPERBOLIC ARC TANGENT OF X 20*24590Szliu * DOUBLE PRECISION (VAX D format 56 bits, IEEE DOUBLE 53 BITS) 21*24590Szliu * CODED IN C BY K.C. NG, 1/8/85; 22*24590Szliu * REVISED BY K.C. NG on 2/7/85, 3/7/85, 8/18/85. 23*24590Szliu * 24*24590Szliu * Required kernel function: 25*24590Szliu * log1p(x) ...return log(1+x) 26*24590Szliu * 27*24590Szliu * Method : 28*24590Szliu * Return 29*24590Szliu * 1 2x x 30*24590Szliu * atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------) 31*24590Szliu * 2 1 - x 1 - x 32*24590Szliu * 33*24590Szliu * Special cases: 34*24590Szliu * atanh(x) is NaN if |x| > 1 with signal; 35*24590Szliu * atanh(NaN) is that NaN with no signal; 36*24590Szliu * atanh(+-1) is +-INF with signal. 37*24590Szliu * 38*24590Szliu * Accuracy: 39*24590Szliu * atanh(x) returns the exact hyperbolic arc tangent of x nearly rounded. 40*24590Szliu * In a test run with 512,000 random arguments on a VAX, the maximum 41*24590Szliu * observed error was 1.87 ulps (units in the last place) at 42*24590Szliu * x= -3.8962076028810414000e-03. 43*24590Szliu */ 44*24590Szliu #ifdef VAX 45*24590Szliu #include <errno.h> 46*24590Szliu #endif 47*24590Szliu 48*24590Szliu double atanh(x) 49*24590Szliu double x; 50*24590Szliu { 51*24590Szliu double copysign(),log1p(),z; 52*24590Szliu z = copysign(0.5,x); 53*24590Szliu x = copysign(x,1.0); 54*24590Szliu #ifdef VAX 55*24590Szliu if (x == 1.0) { 56*24590Szliu extern double infnan(); 57*24590Szliu return(copysign(1.0,z)*infnan(ERANGE)); /* sign(x)*INF */ 58*24590Szliu } 59*24590Szliu #endif 60*24590Szliu x = x/(1.0-x); 61*24590Szliu return( z*log1p(x+x) ); 62*24590Szliu } 63