e_acos.c (7e546392b5fe3a496acff53ac7aadd1c57b2a4cf) e_acos.c (9faa8dc6cc23363b3d8897598cff7d3a40045a46)
1/* @(#)e_acos.c 5.1 93/09/24 */
2/*
3 * ====================================================
4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5 *
6 * Developed at SunPro, a Sun Microsystems, Inc. business.
7 * Permission to use, copy, modify, and distribute this
8 * software is freely granted, provided that this notice
9 * is preserved.
10 * ====================================================
11 */
12
13#ifndef lint
1/* @(#)e_acos.c 5.1 93/09/24 */
2/*
3 * ====================================================
4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5 *
6 * Developed at SunPro, a Sun Microsystems, Inc. business.
7 * Permission to use, copy, modify, and distribute this
8 * software is freely granted, provided that this notice
9 * is preserved.
10 * ====================================================
11 */
12
13#ifndef lint
14static char rcsid[] = "$Id$";
14static char rcsid[] = "$Id: e_acos.c,v 1.5 1997/02/22 15:09:54 peter Exp $";
15#endif
16
17/* __ieee754_acos(x)
18 * Method :
19 * acos(x) = pi/2 - asin(x)
20 * acos(-x) = pi/2 + asin(x)
21 * For |x|<=0.5
22 * acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)

--- 7 unchanged lines hidden (view full) ---

30 * For x<-0.5
31 * acos(x) = pi - 2asin(sqrt((1-|x|)/2))
32 * = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z)
33 *
34 * Special cases:
35 * if x is NaN, return x itself;
36 * if |x|>1, return NaN with invalid signal.
37 *
15#endif
16
17/* __ieee754_acos(x)
18 * Method :
19 * acos(x) = pi/2 - asin(x)
20 * acos(-x) = pi/2 + asin(x)
21 * For |x|<=0.5
22 * acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)

--- 7 unchanged lines hidden (view full) ---

30 * For x<-0.5
31 * acos(x) = pi - 2asin(sqrt((1-|x|)/2))
32 * = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z)
33 *
34 * Special cases:
35 * if x is NaN, return x itself;
36 * if |x|>1, return NaN with invalid signal.
37 *
38 * Function needed: sqrt
38 * Function needed: __ieee754_sqrt
39 */
40
41#include "math.h"
42#include "math_private.h"
43
44#ifdef __STDC__
45static const double
46#else

--- 40 unchanged lines hidden (view full) ---

87 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
88 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
89 r = p/q;
90 return pio2_hi - (x - (pio2_lo-x*r));
91 } else if (hx<0) { /* x < -0.5 */
92 z = (one+x)*0.5;
93 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
94 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
39 */
40
41#include "math.h"
42#include "math_private.h"
43
44#ifdef __STDC__
45static const double
46#else

--- 40 unchanged lines hidden (view full) ---

87 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
88 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
89 r = p/q;
90 return pio2_hi - (x - (pio2_lo-x*r));
91 } else if (hx<0) { /* x < -0.5 */
92 z = (one+x)*0.5;
93 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
94 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
95 s = sqrt(z);
95 s = __ieee754_sqrt(z);
96 r = p/q;
97 w = r*s-pio2_lo;
98 return pi - 2.0*(s+w);
99 } else { /* x > 0.5 */
100 z = (one-x)*0.5;
96 r = p/q;
97 w = r*s-pio2_lo;
98 return pi - 2.0*(s+w);
99 } else { /* x > 0.5 */
100 z = (one-x)*0.5;
101 s = sqrt(z);
101 s = __ieee754_sqrt(z);
102 df = s;
103 SET_LOW_WORD(df,0);
104 c = (z-df*df)/(s+df);
105 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
106 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
107 r = p/q;
108 w = r*s+c;
109 return 2.0*(df+w);
110 }
111}
102 df = s;
103 SET_LOW_WORD(df,0);
104 c = (z-df*df)/(s+df);
105 p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
106 q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
107 r = p/q;
108 w = r*s+c;
109 return 2.0*(df+w);
110 }
111}