e_asin.c (181258) | e_asin.c (218509) |
---|---|
1 2/* @(#)e_asin.c 1.3 95/01/18 */ 3/* 4 * ==================================================== 5 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 6 * 7 * Developed at SunSoft, a Sun Microsystems, Inc. business. 8 * Permission to use, copy, modify, and distribute this 9 * software is freely granted, provided that this notice 10 * is preserved. 11 * ==================================================== 12 */ 13 14#include <sys/cdefs.h> | 1 2/* @(#)e_asin.c 1.3 95/01/18 */ 3/* 4 * ==================================================== 5 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 6 * 7 * Developed at SunSoft, a Sun Microsystems, Inc. business. 8 * Permission to use, copy, modify, and distribute this 9 * software is freely granted, provided that this notice 10 * is preserved. 11 * ==================================================== 12 */ 13 14#include <sys/cdefs.h> |
15__FBSDID("$FreeBSD: head/lib/msun/src/e_asin.c 181258 2008-08-03 17:49:05Z das $"); | 15__FBSDID("$FreeBSD: head/lib/msun/src/e_asin.c 218509 2011-02-10 07:37:50Z das $"); |
16 17/* __ieee754_asin(x) 18 * Method : 19 * Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ... 20 * we approximate asin(x) on [0,0.5] by 21 * asin(x) = x + x*x^2*R(x^2) 22 * where 23 * R(x^2) is a rational approximation of (asin(x)-x)/x^3 --- 53 unchanged lines hidden (view full) --- 77 if(ix>= 0x3ff00000) { /* |x|>= 1 */ 78 u_int32_t lx; 79 GET_LOW_WORD(lx,x); 80 if(((ix-0x3ff00000)|lx)==0) 81 /* asin(1)=+-pi/2 with inexact */ 82 return x*pio2_hi+x*pio2_lo; 83 return (x-x)/(x-x); /* asin(|x|>1) is NaN */ 84 } else if (ix<0x3fe00000) { /* |x|<0.5 */ | 16 17/* __ieee754_asin(x) 18 * Method : 19 * Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ... 20 * we approximate asin(x) on [0,0.5] by 21 * asin(x) = x + x*x^2*R(x^2) 22 * where 23 * R(x^2) is a rational approximation of (asin(x)-x)/x^3 --- 53 unchanged lines hidden (view full) --- 77 if(ix>= 0x3ff00000) { /* |x|>= 1 */ 78 u_int32_t lx; 79 GET_LOW_WORD(lx,x); 80 if(((ix-0x3ff00000)|lx)==0) 81 /* asin(1)=+-pi/2 with inexact */ 82 return x*pio2_hi+x*pio2_lo; 83 return (x-x)/(x-x); /* asin(|x|>1) is NaN */ 84 } else if (ix<0x3fe00000) { /* |x|<0.5 */ |
85 if(ix<0x3e400000) { /* if |x| < 2**-27 */ | 85 if(ix<0x3e500000) { /* if |x| < 2**-26 */ |
86 if(huge+x>one) return x;/* return x with inexact if x!=0*/ 87 } 88 t = x*x; 89 p = t*(pS0+t*(pS1+t*(pS2+t*(pS3+t*(pS4+t*pS5))))); 90 q = one+t*(qS1+t*(qS2+t*(qS3+t*qS4))); 91 w = p/q; 92 return x+x*w; 93 } --- 24 unchanged lines hidden --- | 86 if(huge+x>one) return x;/* return x with inexact if x!=0*/ 87 } 88 t = x*x; 89 p = t*(pS0+t*(pS1+t*(pS2+t*(pS3+t*(pS4+t*pS5))))); 90 q = one+t*(qS1+t*(qS2+t*(qS3+t*qS4))); 91 w = p/q; 92 return x+x*w; 93 } --- 24 unchanged lines hidden --- |