Index: fusl/src/math/sqrtf.c |
diff --git a/fusl/src/math/sqrtf.c b/fusl/src/math/sqrtf.c |
index 28cb4ad3713e6447843631a8d4e686f0f1dfdbb9..329a00a1a06e4759dd12769691d6fd57b4094d7d 100644 |
--- a/fusl/src/math/sqrtf.c |
+++ b/fusl/src/math/sqrtf.c |
@@ -17,68 +17,67 @@ |
static const float tiny = 1.0e-30; |
-float sqrtf(float x) |
-{ |
- float z; |
- int32_t sign = (int)0x80000000; |
- int32_t ix,s,q,m,t,i; |
- uint32_t r; |
+float sqrtf(float x) { |
+ float z; |
+ int32_t sign = (int)0x80000000; |
+ int32_t ix, s, q, m, t, i; |
+ uint32_t r; |
- GET_FLOAT_WORD(ix, x); |
+ GET_FLOAT_WORD(ix, x); |
- /* take care of Inf and NaN */ |
- if ((ix&0x7f800000) == 0x7f800000) |
- return x*x + x; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf, sqrt(-inf)=sNaN */ |
+ /* take care of Inf and NaN */ |
+ if ((ix & 0x7f800000) == 0x7f800000) |
+ return x * x + x; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf, sqrt(-inf)=sNaN */ |
- /* take care of zero */ |
- if (ix <= 0) { |
- if ((ix&~sign) == 0) |
- return x; /* sqrt(+-0) = +-0 */ |
- if (ix < 0) |
- return (x-x)/(x-x); /* sqrt(-ve) = sNaN */ |
- } |
- /* normalize x */ |
- m = ix>>23; |
- if (m == 0) { /* subnormal x */ |
- for (i = 0; (ix&0x00800000) == 0; i++) |
- ix<<=1; |
- m -= i - 1; |
- } |
- m -= 127; /* unbias exponent */ |
- ix = (ix&0x007fffff)|0x00800000; |
- if (m&1) /* odd m, double x to make it even */ |
- ix += ix; |
- m >>= 1; /* m = [m/2] */ |
+ /* take care of zero */ |
+ if (ix <= 0) { |
+ if ((ix & ~sign) == 0) |
+ return x; /* sqrt(+-0) = +-0 */ |
+ if (ix < 0) |
+ return (x - x) / (x - x); /* sqrt(-ve) = sNaN */ |
+ } |
+ /* normalize x */ |
+ m = ix >> 23; |
+ if (m == 0) { /* subnormal x */ |
+ for (i = 0; (ix & 0x00800000) == 0; i++) |
+ ix <<= 1; |
+ m -= i - 1; |
+ } |
+ m -= 127; /* unbias exponent */ |
+ ix = (ix & 0x007fffff) | 0x00800000; |
+ if (m & 1) /* odd m, double x to make it even */ |
+ ix += ix; |
+ m >>= 1; /* m = [m/2] */ |
- /* generate sqrt(x) bit by bit */ |
- ix += ix; |
- q = s = 0; /* q = sqrt(x) */ |
- r = 0x01000000; /* r = moving bit from right to left */ |
+ /* generate sqrt(x) bit by bit */ |
+ ix += ix; |
+ q = s = 0; /* q = sqrt(x) */ |
+ r = 0x01000000; /* r = moving bit from right to left */ |
- while (r != 0) { |
- t = s + r; |
- if (t <= ix) { |
- s = t+r; |
- ix -= t; |
- q += r; |
- } |
- ix += ix; |
- r >>= 1; |
- } |
+ while (r != 0) { |
+ t = s + r; |
+ if (t <= ix) { |
+ s = t + r; |
+ ix -= t; |
+ q += r; |
+ } |
+ ix += ix; |
+ r >>= 1; |
+ } |
- /* use floating add to find out rounding direction */ |
- if (ix != 0) { |
- z = 1.0f - tiny; /* raise inexact flag */ |
- if (z >= 1.0f) { |
- z = 1.0f + tiny; |
- if (z > 1.0f) |
- q += 2; |
- else |
- q += q & 1; |
- } |
- } |
- ix = (q>>1) + 0x3f000000; |
- ix += m << 23; |
- SET_FLOAT_WORD(z, ix); |
- return z; |
+ /* use floating add to find out rounding direction */ |
+ if (ix != 0) { |
+ z = 1.0f - tiny; /* raise inexact flag */ |
+ if (z >= 1.0f) { |
+ z = 1.0f + tiny; |
+ if (z > 1.0f) |
+ q += 2; |
+ else |
+ q += q & 1; |
+ } |
+ } |
+ ix = (q >> 1) + 0x3f000000; |
+ ix += m << 23; |
+ SET_FLOAT_WORD(z, ix); |
+ return z; |
} |