12aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm/** @file 22aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm Compute acos(x) using ieee FP math. 32aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 42aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm Copyright (c) 2010 - 2011, Intel Corporation. All rights reserved.<BR> 52aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm This program and the accompanying materials are licensed and made available under 62aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm the terms and conditions of the BSD License that accompanies this distribution. 72aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm The full text of the license may be found at 82aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm http://opensource.org/licenses/bsd-license. 92aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 102aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm THE PROGRAM IS DISTRIBUTED UNDER THE BSD LICENSE ON AN "AS IS" BASIS, 112aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm WITHOUT WARRANTIES OR REPRESENTATIONS OF ANY KIND, EITHER EXPRESS OR IMPLIED. 122aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 132aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * ==================================================== 142aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 152aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * 162aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Developed at SunPro, a Sun Microsystems, Inc. business. 172aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Permission to use, copy, modify, and distribute this 182aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * software is freely granted, provided that this notice 192aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * is preserved. 202aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * ==================================================== 212aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 222aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm e_acos.c 5.1 93/09/24 232aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm NetBSD: e_acos.c,v 1.12 2002/05/26 22:01:47 wiz Exp 242aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm */ 252aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#if defined(_MSC_VER) /* Handle Microsoft VC++ compiler specifics. */ 262aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm // Keep older compilers quiet about floating-point divide-by-zero 272aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm #pragma warning ( disable : 4723 ) 282aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#endif 292aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 302aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#include <LibConfig.h> 312aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#include <sys/EfiCdefs.h> 322aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 332aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm/* __ieee754_acos(x) 342aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Method : 352aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * acos(x) = pi/2 - asin(x) 362aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * acos(-x) = pi/2 + asin(x) 372aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * For |x|<=0.5 382aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c) 392aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * For x>0.5 402aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * acos(x) = pi/2 - (pi/2 - 2asin(sqrt((1-x)/2))) 412aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * = 2asin(sqrt((1-x)/2)) 422aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * = 2s + 2s*z*R(z) ...z=(1-x)/2, s=sqrt(z) 432aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * = 2f + (2c + 2s*z*R(z)) 442aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * where f=hi part of s, and c = (z-f*f)/(s+f) is the correction term 452aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * for f so that f+c ~ sqrt(z). 462aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * For x<-0.5 472aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * acos(x) = pi - 2asin(sqrt((1-|x|)/2)) 482aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z) 492aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * 502aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Special cases: 512aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * if x is NaN, return x itself; 522aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * if |x|>1, return NaN with invalid signal. 532aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * 542aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm * Function needed: __ieee754_sqrt 552aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm */ 562aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 572aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#include "math.h" 582aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm#include "math_private.h" 592aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 602aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmstatic const double 612aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmone = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */ 622aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */ 632aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */ 642aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpio2_lo = 6.12323399573676603587e-17, /* 0x3C91A626, 0x33145C07 */ 652aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS0 = 1.66666666666666657415e-01, /* 0x3FC55555, 0x55555555 */ 662aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS1 = -3.25565818622400915405e-01, /* 0xBFD4D612, 0x03EB6F7D */ 672aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS2 = 2.01212532134862925881e-01, /* 0x3FC9C155, 0x0E884455 */ 682aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS3 = -4.00555345006794114027e-02, /* 0xBFA48228, 0xB5688F3B */ 692aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS4 = 7.91534994289814532176e-04, /* 0x3F49EFE0, 0x7501B288 */ 702aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmpS5 = 3.47933107596021167570e-05, /* 0x3F023DE1, 0x0DFDF709 */ 712aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmqS1 = -2.40339491173441421878e+00, /* 0xC0033A27, 0x1C8A2D4B */ 722aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmqS2 = 2.02094576023350569471e+00, /* 0x40002AE5, 0x9C598AC8 */ 732aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmqS3 = -6.88283971605453293030e-01, /* 0xBFE6066C, 0x1B8D0159 */ 742aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmqS4 = 7.70381505559019352791e-02; /* 0x3FB3B8C5, 0xB12E9282 */ 752aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 762aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylmdouble 772aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm__ieee754_acos(double x) 782aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm{ 792aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm double z,p,q,r,w,s,c,df; 802aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm int32_t hx,ix; 812aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm GET_HIGH_WORD(hx,x); 822aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm ix = hx&0x7fffffff; 832aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm if(ix>=0x3ff00000) { /* |x| >= 1 */ 842aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm u_int32_t lx; 852aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm 862aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm GET_LOW_WORD(lx,x); 872aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm if(((ix-0x3ff00000)|lx)==0) { /* |x|==1 */ 882aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm if(hx>0) return 0.0; /* acos(1) = 0 */ 892aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm else return pi+2.0*pio2_lo; /* acos(-1)= pi */ 902aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm } 912aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm return (x-x)/(x-x); /* acos(|x|>1) is NaN */ 922aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm } 932aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm if(ix<0x3fe00000) { /* |x| < 0.5 */ 942aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm if(ix<=0x3c600000) return pio2_hi+pio2_lo; /*if|x|<2**-57*/ 952aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm z = x*x; 962aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5))))); 972aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4))); 982aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm r = p/q; 992aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm return pio2_hi - (x - (pio2_lo-x*r)); 1002aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm } 1012aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm else if (hx<0) { /* x < -0.5 */ 1022aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm z = (one+x)*0.5; 1032aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5))))); 1042aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4))); 1052aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm s = __ieee754_sqrt(z); 1062aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm r = p/q; 1072aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm w = r*s-pio2_lo; 1082aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm return pi - 2.0*(s+w); 1092aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm } 1102aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm else { /* x > 0.5 */ 1112aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm z = (one-x)*0.5; 1122aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm s = __ieee754_sqrt(z); 1132aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm df = s; 1142aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm SET_LOW_WORD(df,0); 1152aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm c = (z-df*df)/(s+df); 1162aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5))))); 1172aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4))); 1182aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm r = p/q; 1192aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm w = r*s+c; 1202aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm return 2.0*(df+w); 1212aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm } 1222aa62f2bc9a9654687b377d9ca8a8c2c860a3852darylm} 123