Deleted Added
sdiff udiff text old ( 176467 ) new ( 176476 )
full compact
e_rem_pio2l.h (176467) e_rem_pio2l.h (176476)
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 *
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 * Optimized by Bruce D. Evans.
13 */
14
15#include <sys/cdefs.h>
14 */
15
16#include <sys/cdefs.h>
16__FBSDID("$FreeBSD: head/lib/msun/src/e_rem_pio2.c 176467 2008-02-22 18:43:23Z bde $");
17__FBSDID("$FreeBSD: head/lib/msun/src/e_rem_pio2.c 176476 2008-02-23 12:53:21Z 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>

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

122 y[0] = z + 4*pio2_1t;
123 y[1] = (z-y[0])+4*pio2_1t;
124 return -4;
125 }
126 }
127 }
128 if(ix<=0x413921fb) { /* |x| ~<= 2^19*(pi/2), medium size */
129medium:
18
19/* __ieee754_rem_pio2(x,y)
20 *
21 * return the remainder of x rem pi/2 in y[0]+y[1]
22 * use __kernel_rem_pio2()
23 */
24
25#include <float.h>

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

123 y[0] = z + 4*pio2_1t;
124 y[1] = (z-y[0])+4*pio2_1t;
125 return -4;
126 }
127 }
128 }
129 if(ix<=0x413921fb) { /* |x| ~<= 2^19*(pi/2), medium size */
130medium:
130 t = fabs(x);
131 /* Use a specialized rint() to get fn. Assume round-to-nearest. */
131 /* Use a specialized rint() to get fn. Assume round-to-nearest. */
132 STRICT_ASSIGN(double,fn,t*invpio2+0x1.8p52);
132 STRICT_ASSIGN(double,fn,x*invpio2+0x1.8p52);
133 fn = fn-0x1.8p52;
134#ifdef HAVE_EFFICIENT_IRINT
135 n = irint(fn);
136#else
137 n = (int32_t)fn;
138#endif
133 fn = fn-0x1.8p52;
134#ifdef HAVE_EFFICIENT_IRINT
135 n = irint(fn);
136#else
137 n = (int32_t)fn;
138#endif
139 r = t-fn*pio2_1;
139 r = x-fn*pio2_1;
140 w = fn*pio2_1t; /* 1st round good to 85 bit */
141 {
142 u_int32_t high;
143 j = ix>>20;
144 y[0] = r-w;
145 GET_HIGH_WORD(high,y[0]);
146 i = j-((high>>20)&0x7ff);
147 if(i>16) { /* 2nd iteration needed, good to 118 */

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

157 w = fn*pio2_3;
158 r = t-w;
159 w = fn*pio2_3t-((t-r)-w);
160 y[0] = r-w;
161 }
162 }
163 }
164 y[1] = (r-y[0])-w;
140 w = fn*pio2_1t; /* 1st round good to 85 bit */
141 {
142 u_int32_t high;
143 j = ix>>20;
144 y[0] = r-w;
145 GET_HIGH_WORD(high,y[0]);
146 i = j-((high>>20)&0x7ff);
147 if(i>16) { /* 2nd iteration needed, good to 118 */

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

157 w = fn*pio2_3;
158 r = t-w;
159 w = fn*pio2_3t-((t-r)-w);
160 y[0] = r-w;
161 }
162 }
163 }
164 y[1] = (r-y[0])-w;
165 if(hx<0) {y[0] = -y[0]; y[1] = -y[1]; return -n;}
166 else return n;
165 return n;
167 }
168 /*
169 * all other (large) arguments
170 */
171 if(ix>=0x7ff00000) { /* x is inf or NaN */
172 y[0]=y[1]=x-x; return 0;
173 }
174 /* set z = scalbn(|x|,ilogb(x)-23) */

--- 15 unchanged lines hidden --- 		166 }
167 /*
168 * all other (large) arguments
169 */
170 if(ix>=0x7ff00000) { /* x is inf or NaN */
171 y[0]=y[1]=x-x; return 0;
172 }
173 /* set z = scalbn(|x|,ilogb(x)-23) */

--- 15 unchanged lines hidden ---