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 __asinl = asinl 31 32 /* 33 * asinl(x) = atan2l(x,sqrt(1-x*x)); 34 * 35 * For better accuracy, 1-x*x is computed as follows 36 * 1-x*x if x < 0.5, 37 * 2*(1-|x|)-(1-|x|)*(1-|x|) if x >= 0.5. 38 * 39 * Special cases: 40 * if x is NaN, return x itself; 41 * if |x|>1, return NaN with invalid signal. 42 */ 43 44 #include "libm.h" 45 46 static const long double zero = 0.0L, small = 1.0e-20L, half = 0.5L, one = 1.0L; 47 #ifndef lint 48 static const long double big = 1.0e+20L; 49 #endif 50 51 long double 52 asinl(long double x) { 53 long double t, w; 54 volatile long double dummy; 55 56 w = fabsl(x); 57 if (isnanl(x)) 58 return (x + x); 59 else if (w <= half) { 60 if (w < small) { 61 #ifndef lint 62 dummy = w + big; 63 /* inexact if w != 0 */ 64 #endif 65 return (x); 66 } else 67 return (atanl(x / sqrtl(one - x * x))); 68 } else if (w < one) { 69 t = one - w; 70 w = t + t; 71 return (atanl(x / sqrtl(w - t * t))); 72 } else if (w == one) 73 return (atan2l(x, zero)); /* asin(+-1) = +- PI/2 */ 74 else 75 return (zero / zero); /* |x| > 1: invalid */ 76 } 77