1 /* e_sqrtf.c -- float version of e_sqrt.c.
5 * ====================================================
6 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
8 * Developed at SunPro, a Sun Microsystems, Inc. business.
9 * Permission to use, copy, modify, and distribute this
10 * software is freely granted, provided that this notice
12 * ====================================================
16 #include <math_private.h>
17 #include <libm-alias-finite.h>
18 #include <math-use-builtins.h>
21 __ieee754_sqrtf(float x
)
24 return __builtin_sqrtf (x
);
26 /* Use generic implementation. */
28 int32_t sign
= (int)0x80000000;
34 /* take care of Inf and NaN */
35 if((ix
&0x7f800000)==0x7f800000) {
36 return x
*x
+x
; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf
39 /* take care of zero */
41 if((ix
&(~sign
))==0) return x
;/* sqrt(+-0) = +-0 */
43 return (x
-x
)/(x
-x
); /* sqrt(-ve) = sNaN */
47 if(m
==0) { /* subnormal x */
48 for(i
=0;(ix
&0x00800000)==0;i
++) ix
<<=1;
51 m
-= 127; /* unbias exponent */
52 ix
= (ix
&0x007fffff)|0x00800000;
53 if(m
&1) /* odd m, double x to make it even */
55 m
>>= 1; /* m = [m/2] */
57 /* generate sqrt(x) bit by bit */
59 q
= s
= 0; /* q = sqrt(x) */
60 r
= 0x01000000; /* r = moving bit from right to left */
73 /* use floating add to find out rounding direction */
75 z
= 0x1p
0 - 0x1.4484cp
-100; /* trigger inexact flag. */
77 z
= 0x1p
0 + 0x1.4484cp
-100;
84 ix
= (q
>>1)+0x3f000000;
88 #endif /* ! USE_SQRTF_BUILTIN */
90 #ifndef __ieee754_sqrtf
91 libm_alias_finite (__ieee754_sqrtf
, __sqrtf
)