4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync/** @file
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync Compute acos(x) using ieee FP math.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync Copyright (c) 2010 - 2011, Intel Corporation. All rights reserved.<BR>
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync This program and the accompanying materials are licensed and made available under
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync the terms and conditions of the BSD License that accompanies this distribution.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync The full text of the license may be found at
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync http://opensource.org/licenses/bsd-license.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync THE PROGRAM IS DISTRIBUTED UNDER THE BSD LICENSE ON AN "AS IS" BASIS,
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync WITHOUT WARRANTIES OR REPRESENTATIONS OF ANY KIND, EITHER EXPRESS OR IMPLIED.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * ====================================================
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync *
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Developed at SunPro, a Sun Microsystems, Inc. business.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Permission to use, copy, modify, and distribute this
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * software is freely granted, provided that this notice
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * is preserved.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * ====================================================
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync e_acos.c 5.1 93/09/24
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync NetBSD: e_acos.c,v 1.12 2002/05/26 22:01:47 wiz Exp
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#if defined(_MSC_VER) /* Handle Microsoft VC++ compiler specifics. */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync // Keep older compilers quiet about floating-point divide-by-zero
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync #pragma warning ( disable : 4723 )
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#endif
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#include <LibConfig.h>
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#include <sys/EfiCdefs.h>
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync/* __ieee754_acos(x)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Method :
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * acos(x) = pi/2 - asin(x)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * acos(-x) = pi/2 + asin(x)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * For |x|<=0.5
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * For x>0.5
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * acos(x) = pi/2 - (pi/2 - 2asin(sqrt((1-x)/2)))
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * = 2asin(sqrt((1-x)/2))
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * = 2s + 2s*z*R(z) ...z=(1-x)/2, s=sqrt(z)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * = 2f + (2c + 2s*z*R(z))
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * where f=hi part of s, and c = (z-f*f)/(s+f) is the correction term
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * for f so that f+c ~ sqrt(z).
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * For x<-0.5
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * acos(x) = pi - 2asin(sqrt((1-|x|)/2))
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync *
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Special cases:
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * if x is NaN, return x itself;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * if |x|>1, return NaN with invalid signal.
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync *
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync * Function needed: __ieee754_sqrt
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#include "math.h"
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync#include "math_private.h"
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncstatic const double
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncone = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpio2_lo = 6.12323399573676603587e-17, /* 0x3C91A626, 0x33145C07 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS0 = 1.66666666666666657415e-01, /* 0x3FC55555, 0x55555555 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS1 = -3.25565818622400915405e-01, /* 0xBFD4D612, 0x03EB6F7D */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS2 = 2.01212532134862925881e-01, /* 0x3FC9C155, 0x0E884455 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS3 = -4.00555345006794114027e-02, /* 0xBFA48228, 0xB5688F3B */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS4 = 7.91534994289814532176e-04, /* 0x3F49EFE0, 0x7501B288 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncpS5 = 3.47933107596021167570e-05, /* 0x3F023DE1, 0x0DFDF709 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncqS1 = -2.40339491173441421878e+00, /* 0xC0033A27, 0x1C8A2D4B */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncqS2 = 2.02094576023350569471e+00, /* 0x40002AE5, 0x9C598AC8 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncqS3 = -6.88283971605453293030e-01, /* 0xBFE6066C, 0x1B8D0159 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncqS4 = 7.70381505559019352791e-02; /* 0x3FB3B8C5, 0xB12E9282 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsyncdouble
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync__ieee754_acos(double x)
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync{
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync double z,p,q,r,w,s,c,df;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync int32_t hx,ix;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync GET_HIGH_WORD(hx,x);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync ix = hx&0x7fffffff;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync if(ix>=0x3ff00000) { /* |x| >= 1 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync u_int32_t lx;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync GET_LOW_WORD(lx,x);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync if(((ix-0x3ff00000)|lx)==0) { /* |x|==1 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync if(hx>0) return 0.0; /* acos(1) = 0 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync else return pi+2.0*pio2_lo; /* acos(-1)= pi */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync }
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync return (x-x)/(x-x); /* acos(|x|>1) is NaN */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync }
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync if(ix<0x3fe00000) { /* |x| < 0.5 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync if(ix<=0x3c600000) return pio2_hi+pio2_lo; /*if|x|<2**-57*/
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync z = x*x;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync r = p/q;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync return pio2_hi - (x - (pio2_lo-x*r));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync }
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync else if (hx<0) { /* x < -0.5 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync z = (one+x)*0.5;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync s = __ieee754_sqrt(z);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync r = p/q;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync w = r*s-pio2_lo;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync return pi - 2.0*(s+w);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync }
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync else { /* x > 0.5 */
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync z = (one-x)*0.5;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync s = __ieee754_sqrt(z);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync df = s;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync SET_LOW_WORD(df,0);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync c = (z-df*df)/(s+df);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync r = p/q;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync w = r*s+c;
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync return 2.0*(df+w);
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync }
4fd606d1f5abe38e1f42c38de1d2e895166bd0f4vboxsync}