Print this page
5262 libm needs to be carefully unifdef'd
5268 libm doesn't need to hide symbols which are already local
@@ -28,82 +28,13 @@
#pragma weak sqrtf = __sqrtf
#include "libm.h"
-#ifdef __INLINE
extern float __inline_sqrtf(float);
float
sqrtf(float x) {
return (__inline_sqrtf(x));
}
-#else /* defined(__INLINE) */
-
-static const float huge = 1.0e35F, tiny = 1.0e-35F, zero = 0.0f;
-
-float
-sqrtf(float x) {
- float dz, w;
- int *pw = (int *)&w;
- int ix, j, r, q, m, n, s, t;
-
- w = x;
- ix = pw[0];
- if (ix <= 0) {
- /* x is <= 0 or nan */
- j = ix & 0x7fffffff;
- if (j == 0)
- return (w);
- return ((w * zero) / zero);
- }
-
- if ((ix & 0x7f800000) == 0x7f800000) {
- /* x is +inf or nan */
- return (w * w);
- }
-
- m = ir_ilogb_(&w);
- n = -m;
- w = r_scalbn_(&w, (int *)&n);
- ix = (pw[0] & 0x007fffff) | 0x00800000;
- n = m / 2;
- if ((n + n) != m) {
- ix = ix + ix;
- m -= 1;
- n = m / 2;
- }
-
- /* generate sqrt(x) bit by bit */
- ix <<= 1;
- q = s = 0;
- r = 0x01000000;
- for (j = 1; j <= 25; j++) {
- t = s + r;
- if (t <= ix) {
- s = t + r;
- ix -= t;
- q += r;
- }
- ix <<= 1;
- r >>= 1;
- }
- if (ix == 0)
- goto done;
-
- /* raise inexact and determine the ambient rounding mode */
- dz = huge - tiny;
- if (dz < huge)
- goto done;
- dz = huge + tiny;
- if (dz > huge)
- q += 1;
- q += (q & 1);
-
-done:
- pw[0] = (q >> 1) + 0x3f000000;
- return (r_scalbn_(&w, (int *)&n));
-}
-
-#endif /* defined(__INLINE) */