12116Sjkh/* @(#)s_asinh.c 5.1 93/09/24 */
22116Sjkh/*
32116Sjkh * ====================================================
42116Sjkh * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
52116Sjkh *
62116Sjkh * Developed at SunPro, a Sun Microsystems, Inc. business.
72116Sjkh * Permission to use, copy, modify, and distribute this
88870Srgrimes * software is freely granted, provided that this notice
92116Sjkh * is preserved.
102116Sjkh * ====================================================
112116Sjkh */
122116Sjkh
13176451Sdas#include <sys/cdefs.h>
14176451Sdas__FBSDID("$FreeBSD: releng/11.0/lib/msun/src/s_asinh.c 251599 2013-06-10 06:04:58Z das $");
152116Sjkh
162116Sjkh/* asinh(x)
172116Sjkh * Method :
188870Srgrimes *	Based on
192116Sjkh *		asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
202116Sjkh *	we have
212116Sjkh *	asinh(x) := x  if  1+x*x=1,
222116Sjkh *		 := sign(x)*(log(x)+ln2)) for large |x|, else
232116Sjkh *		 := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
248870Srgrimes *		 := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
252116Sjkh */
262116Sjkh
27251599Sdas#include <float.h>
28251599Sdas
292116Sjkh#include "math.h"
302116Sjkh#include "math_private.h"
312116Sjkh
328870Srgrimesstatic const double
332116Sjkhone =  1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
342116Sjkhln2 =  6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
358870Srgrimeshuge=  1.00000000000000000000e+300;
362116Sjkh
3797413Salfreddouble
3897413Salfredasinh(double x)
398870Srgrimes{
402116Sjkh	double t,w;
412116Sjkh	int32_t hx,ix;
422116Sjkh	GET_HIGH_WORD(hx,x);
432116Sjkh	ix = hx&0x7fffffff;
442116Sjkh	if(ix>=0x7ff00000) return x+x;	/* x is inf or NaN */
452116Sjkh	if(ix< 0x3e300000) {	/* |x|<2**-28 */
462116Sjkh	    if(huge+x>one) return x;	/* return x inexact except 0 */
478870Srgrimes	}
482116Sjkh	if(ix>0x41b00000) {	/* |x| > 2**28 */
492116Sjkh	    w = __ieee754_log(fabs(x))+ln2;
502116Sjkh	} else if (ix>0x40000000) {	/* 2**28 > |x| > 2.0 */
512116Sjkh	    t = fabs(x);
5223579Sbde	    w = __ieee754_log(2.0*t+one/(__ieee754_sqrt(x*x+one)+t));
532116Sjkh	} else {		/* 2.0 > |x| > 2**-28 */
542116Sjkh	    t = x*x;
5523579Sbde	    w =log1p(fabs(x)+t/(one+__ieee754_sqrt(one+t)));
562116Sjkh	}
572116Sjkh	if(hx>0) return w; else return -w;
582116Sjkh}
59251599Sdas
60251599Sdas#if LDBL_MANT_DIG == 53
61251599Sdas__weak_reference(asinh, asinhl);
62251599Sdas#endif
63