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 acosh = __acosh 31 32 /* INDENT OFF */ 33 /* 34 * acosh(x) 35 * Method : 36 * Based on 37 * acosh(x) = log [ x + sqrt(x*x-1) ] 38 * we have 39 * acosh(x) := log(x)+ln2, if x is large; else 40 * acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x > 2; else 41 * acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t = x-1. 42 * 43 * Special cases: 44 * acosh(x) is NaN with signal if x < 1. 45 * acosh(NaN) is NaN without signal. 46 */ 47 /* INDENT ON */ 48 49 #include "libm_synonyms.h" /* __acosh, __log, __log1p */ 50 #include "libm_protos.h" /* _SVID_libm_error */ 51 #include "libm_macros.h" 52 #include <math.h> 53 54 static const double 55 one = 1.0, 56 ln2 = 6.93147180559945286227e-01; /* 3FE62E42, FEFA39EF */ 57 58 double 59 acosh(double x) { 60 double t; 61 int hx; 62 63 hx = ((int *) &x)[HIWORD]; 64 if (hx < 0x3ff00000) { /* x < 1 */ 65 if (isnan(x)) 66 #if defined(FPADD_TRAPS_INCOMPLETE_ON_NAN) 67 return (hx >= 0xfff80000 ? x : (x - x) / (x - x)); 68 /* assumes sparc-like QNaN */ 69 #else 70 return (x - x) / (x - x); 71 #endif 72 else 73 return (_SVID_libm_err(x, x, 29)); 74 } else if (hx >= 0x41b00000) { 75 /* x > 2**28 */ 76 if (hx >= 0x7ff00000) { /* x is inf of NaN */ 77 #if defined(FPADD_TRAPS_INCOMPLETE_ON_NAN) 78 return (hx >= 0x7ff80000 ? x : x + x); 79 /* assumes sparc-like QNaN */ 80 #else 81 return (x + x); 82 #endif 83 } else /* acosh(huge)=log(2x) */ 84 return (log(x) + ln2); 85 } else if (((hx - 0x3ff00000) | ((int *) &x)[LOWORD]) == 0) { 86 return (0.0); /* acosh(1) = 0 */ 87 } else if (hx > 0x40000000) { 88 /* 2**28 > x > 2 */ 89 t = x * x; 90 return (log(2.0 * x - one / (x + sqrt(t - one)))); 91 } else { 92 /* 1 < x < 2 */ 93 t = x - one; 94 return (log1p(t + sqrt(2.0 * t + t * t))); 95 } 96 } 97