| OLD | NEW |
| 1 #include "libm.h" | 1 #include "libm.h" |
| 2 | 2 |
| 3 #if FLT_EVAL_METHOD==0 || FLT_EVAL_METHOD==1 | 3 #if FLT_EVAL_METHOD == 0 || FLT_EVAL_METHOD == 1 |
| 4 #define EPS DBL_EPSILON | 4 #define EPS DBL_EPSILON |
| 5 #elif FLT_EVAL_METHOD==2 | 5 #elif FLT_EVAL_METHOD == 2 |
| 6 #define EPS LDBL_EPSILON | 6 #define EPS LDBL_EPSILON |
| 7 #endif | 7 #endif |
| 8 static const double_t toint = 1/EPS; | 8 static const double_t toint = 1 / EPS; |
| 9 | 9 |
| 10 double floor(double x) | 10 double floor(double x) { |
| 11 { | 11 union { |
| 12 » union {double f; uint64_t i;} u = {x}; | 12 double f; |
| 13 » int e = u.i >> 52 & 0x7ff; | 13 uint64_t i; |
| 14 » double_t y; | 14 } u = {x}; |
| 15 int e = u.i >> 52 & 0x7ff; |
| 16 double_t y; |
| 15 | 17 |
| 16 » if (e >= 0x3ff+52 || x == 0) | 18 if (e >= 0x3ff + 52 || x == 0) |
| 17 » » return x; | 19 return x; |
| 18 » /* y = int(x) - x, where int(x) is an integer neighbor of x */ | 20 /* y = int(x) - x, where int(x) is an integer neighbor of x */ |
| 19 » if (u.i >> 63) | 21 if (u.i >> 63) |
| 20 » » y = x - toint + toint - x; | 22 y = x - toint + toint - x; |
| 21 » else | 23 else |
| 22 » » y = x + toint - toint - x; | 24 y = x + toint - toint - x; |
| 23 » /* special case because of non-nearest rounding modes */ | 25 /* special case because of non-nearest rounding modes */ |
| 24 » if (e <= 0x3ff-1) { | 26 if (e <= 0x3ff - 1) { |
| 25 » » FORCE_EVAL(y); | 27 FORCE_EVAL(y); |
| 26 » » return u.i >> 63 ? -1 : 0; | 28 return u.i >> 63 ? -1 : 0; |
| 27 » } | 29 } |
| 28 » if (y > 0) | 30 if (y > 0) |
| 29 » » return x + y - 1; | 31 return x + y - 1; |
| 30 » return x + y; | 32 return x + y; |
| 31 } | 33 } |
| OLD | NEW |