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 |