| 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 round(double x) | 10 double round(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) | 18 if (e >= 0x3ff + 52) |
| 17 » » return x; | 19 return x; |
| 18 » if (u.i >> 63) | 20 if (u.i >> 63) |
| 19 » » x = -x; | 21 x = -x; |
| 20 » if (e < 0x3ff-1) { | 22 if (e < 0x3ff - 1) { |
| 21 » » /* raise inexact if x!=0 */ | 23 /* raise inexact if x!=0 */ |
| 22 » » FORCE_EVAL(x + toint); | 24 FORCE_EVAL(x + toint); |
| 23 » » return 0*u.f; | 25 return 0 * u.f; |
| 24 » } | 26 } |
| 25 » y = x + toint - toint - x; | 27 y = x + toint - toint - x; |
| 26 » if (y > 0.5) | 28 if (y > 0.5) |
| 27 » » y = y + x - 1; | 29 y = y + x - 1; |
| 28 » else if (y <= -0.5) | 30 else if (y <= -0.5) |
| 29 » » y = y + x + 1; | 31 y = y + x + 1; |
| 30 » else | 32 else |
| 31 » » y = y + x; | 33 y = y + x; |
| 32 » if (u.i >> 63) | 34 if (u.i >> 63) |
| 33 » » y = -y; | 35 y = -y; |
| 34 » return y; | 36 return y; |
| 35 } | 37 } |
| OLD | NEW |