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 __asinh = asinh 31 32 /* INDENT OFF */ 33 /* 34 * asinh(x) 35 * Method : 36 * Based on 37 * asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ] 38 * we have 39 * asinh(x) := x if 1+x*x == 1, 40 * := sign(x)*(log(x)+ln2)) for large |x|, else 41 * := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x| > 2, else 42 * := sign(x)*log1p(|x|+x^2/(1+sqrt(1+x^2))) 43 */ 44 /* INDENT ON */ 45 46 #include "libm_macros.h" 47 #include <math.h> 48 49 static const double xxx[] = { 50 /* one */ 1.00000000000000000000e+00, /* 3FF00000, 00000000 */ 51 /* ln2 */ 6.93147180559945286227e-01, /* 3FE62E42, FEFA39EF */ 52 /* huge */ 1.00000000000000000000e+300 53 }; 54 #define one xxx[0] 55 #define ln2 xxx[1] 56 #define huge xxx[2] 57 58 double 59 asinh(double x) { 60 double t, w; 61 int hx, ix; 62 63 hx = ((int *) &x)[HIWORD]; 64 ix = hx & 0x7fffffff; 65 if (ix >= 0x7ff00000) 66 #if defined(FPADD_TRAPS_INCOMPLETE_ON_NAN) 67 return (ix >= 0x7ff80000 ? x : x + x); 68 /* assumes sparc-like QNaN */ 69 #else 70 return (x + x); /* x is inf or NaN */ 71 #endif 72 if (ix < 0x3e300000) { /* |x|<2**-28 */ 73 if (huge + x > one) 74 return (x); /* return x inexact except 0 */ 75 } 76 if (ix > 0x41b00000) { /* |x| > 2**28 */ 77 w = log(fabs(x)) + ln2; 78 } else if (ix > 0x40000000) { 79 /* 2**28 > |x| > 2.0 */ 80 t = fabs(x); 81 w = log(2.0 * t + one / (sqrt(x * x + one) + t)); 82 } else { 83 /* 2.0 > |x| > 2**-28 */ 84 t = x * x; 85 w = log1p(fabs(x) + t / (one + sqrt(one + t))); 86 } 87 return (hx > 0 ? w : -w); 88 } 89