Deleted Added
sdiff udiff text old ( 22993 ) new ( 23579 )
full compact
e_hypot.c (22993) e_hypot.c (23579)
1/* @(#)e_hypot.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_hypot.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_hypot.c,v 1.4 1997/02/22 15:10:12 peter Exp $";
15#endif
16
17/* __ieee754_hypot(x,y)
18 *
19 * Method :
20 * If (assume round-to-nearest) z=x*x+y*y
21 * has error less than sqrt(2)/2 ulp, than
22 * sqrt(z) has error less than 1 ulp (exercise).

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

102 }
103 }
104 /* medium size a and b */
105 w = a-b;
106 if (w>b) {
107 t1 = 0;
108 SET_HIGH_WORD(t1,ha);
109 t2 = a-t1;
15#endif
16
17/* __ieee754_hypot(x,y)
18 *
19 * Method :
20 * If (assume round-to-nearest) z=x*x+y*y
21 * has error less than sqrt(2)/2 ulp, than
22 * sqrt(z) has error less than 1 ulp (exercise).

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

102 }
103 }
104 /* medium size a and b */
105 w = a-b;
106 if (w>b) {
107 t1 = 0;
108 SET_HIGH_WORD(t1,ha);
109 t2 = a-t1;
110 w = sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
110 w = __ieee754_sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
111 } else {
112 a = a+a;
113 y1 = 0;
114 SET_HIGH_WORD(y1,hb);
115 y2 = b - y1;
116 t1 = 0;
117 SET_HIGH_WORD(t1,ha+0x00100000);
118 t2 = a - t1;
111 } else {
112 a = a+a;
113 y1 = 0;
114 SET_HIGH_WORD(y1,hb);
115 y2 = b - y1;
116 t1 = 0;
117 SET_HIGH_WORD(t1,ha+0x00100000);
118 t2 = a - t1;
119 w = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
119 w = __ieee754_sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
120 }
121 if(k!=0) {
122 u_int32_t high;
123 t1 = 1.0;
124 GET_HIGH_WORD(high,t1);
125 SET_HIGH_WORD(t1,high+(k<<20));
126 return t1*w;
127 } else return w;
128}
120 }
121 if(k!=0) {
122 u_int32_t high;
123 t1 = 1.0;
124 GET_HIGH_WORD(high,t1);
125 SET_HIGH_WORD(t1,high+(k<<20));
126 return t1*w;
127 } else return w;
128}