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