Deleted Added
sdiff udiff text old ( 176465 ) new ( 176466 )
full compact
e_rem_pio2.c (176465) e_rem_pio2.c (176466)
1
2/* @(#)e_rem_pio2.c 1.4 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
15#include <sys/cdefs.h>
1
2/* @(#)e_rem_pio2.c 1.4 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
15#include <sys/cdefs.h>
16__FBSDID("$FreeBSD: head/lib/msun/src/e_rem_pio2.c 176465 2008-02-22 15:55:14Z bde $");
16__FBSDID("$FreeBSD: head/lib/msun/src/e_rem_pio2.c 176466 2008-02-22 17:26:24Z bde $");
17
18/* __ieee754_rem_pio2(x,y)
19 *
20 * return the remainder of x rem pi/2 in y[0]+y[1]
21 * use __kernel_rem_pio2()
22 */
23
24#include <float.h>
25
26#include "math.h"
27#include "math_private.h"
28
17
18/* __ieee754_rem_pio2(x,y)
19 *
20 * return the remainder of x rem pi/2 in y[0]+y[1]
21 * use __kernel_rem_pio2()
22 */
23
24#include <float.h>
25
26#include "math.h"
27#include "math_private.h"
28
29static const int32_t npio2_hw[] = {
300x3FF921FB, 0x400921FB, 0x4012D97C, 0x401921FB, 0x401F6A7A, 0x4022D97C,
310x4025FDBB, 0x402921FB, 0x402C463A, 0x402F6A7A, 0x4031475C, 0x4032D97C,
320x40346B9C, 0x4035FDBB, 0x40378FDB, 0x403921FB, 0x403AB41B, 0x403C463A,
330x403DD85A, 0x403F6A7A, 0x40407E4C, 0x4041475C, 0x4042106C, 0x4042D97C,
340x4043A28C, 0x40446B9C, 0x404534AC, 0x4045FDBB, 0x4046C6CB, 0x40478FDB,
350x404858EB, 0x404921FB,
36};
37
38/*
39 * invpio2: 53 bits of 2/pi
40 * pio2_1: first 33 bit of pi/2
41 * pio2_1t: pi/2 - pio2_1
42 * pio2_2: second 33 bit of pi/2
43 * pio2_2t: pi/2 - (pio2_1+pio2_2)
44 * pio2_3: third 33 bit of pi/2
45 * pio2_3t: pi/2 - (pio2_1+pio2_2+pio2_3)

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

143 fn = fn-0x1.8p52;
144 n = irint(fn);
145#else
146 n = (int32_t) (t*invpio2+half);
147 fn = (double)n;
148#endif
149 r = t-fn*pio2_1;
150 w = fn*pio2_1t; /* 1st round good to 85 bit */
29/*
30 * invpio2: 53 bits of 2/pi
31 * pio2_1: first 33 bit of pi/2
32 * pio2_1t: pi/2 - pio2_1
33 * pio2_2: second 33 bit of pi/2
34 * pio2_2t: pi/2 - (pio2_1+pio2_2)
35 * pio2_3: third 33 bit of pi/2
36 * pio2_3t: pi/2 - (pio2_1+pio2_2+pio2_3)

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

134 fn = fn-0x1.8p52;
135 n = irint(fn);
136#else
137 n = (int32_t) (t*invpio2+half);
138 fn = (double)n;
139#endif
140 r = t-fn*pio2_1;
141 w = fn*pio2_1t; /* 1st round good to 85 bit */
151 if(n<32&&ix!=npio2_hw[n-1]) {
152 y[0] = r-w; /* quick check no cancellation */
153 } else {
142 {
154 u_int32_t high;
155 j = ix>>20;
156 y[0] = r-w;
157 GET_HIGH_WORD(high,y[0]);
158 i = j-((high>>20)&0x7ff);
159 if(i>16) { /* 2nd iteration needed, good to 118 */
160 t = r;
161 w = fn*pio2_2;

--- 40 unchanged lines hidden --- 		143 u_int32_t high;
144 j = ix>>20;
145 y[0] = r-w;
146 GET_HIGH_WORD(high,y[0]);
147 i = j-((high>>20)&0x7ff);
148 if(i>16) { /* 2nd iteration needed, good to 118 */
149 t = r;
150 w = fn*pio2_2;

--- 40 unchanged lines hidden ---