[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
master 5e229f8 1/2: Fix rounding errors in time conversion
From: |
Paul Eggert |
Subject: |
master 5e229f8 1/2: Fix rounding errors in time conversion |
Date: |
Tue, 3 Mar 2020 13:20:38 -0500 (EST) |
branch: master
commit 5e229f88f946c9c246966defeb629965f65d9549
Author: Paul Eggert <address@hidden>
Commit: Paul Eggert <address@hidden>
Fix rounding errors in time conversion
* src/timefns.c (frac_to_double): Pass FLT_RADIX to mpz_sizeinbase
instead of doing the radix calculation ourselves, not always
correctly. Fix off-by-one error in scale, which caused
double-rounding.
(decode_time_components): Use frac_to_double (via decode_ticks_hz)
to fix double-rounding error that can occur even though
intermediate results are long double.
* test/src/timefns-tests.el (float-time-precision):
Test the above fixes.
---
src/timefns.c | 79 +++++++++++++++++++++--------------------------
test/src/timefns-tests.el | 3 ++
2 files changed, 38 insertions(+), 44 deletions(-)
diff --git a/src/timefns.c b/src/timefns.c
index b301c58..b96a6cb 100644
--- a/src/timefns.c
+++ b/src/timefns.c
@@ -576,18 +576,14 @@ frac_to_double (Lisp_Object numerator, Lisp_Object
denominator)
mpz_t const *n = bignum_integer (&mpz[0], numerator);
mpz_t const *d = bignum_integer (&mpz[1], denominator);
- ptrdiff_t nbits = mpz_sizeinbase (*n, 2);
- ptrdiff_t dbits = mpz_sizeinbase (*d, 2);
- eassume (0 < nbits);
- eassume (0 < dbits);
- ptrdiff_t ndig = (nbits + LOG2_FLT_RADIX - 1) / LOG2_FLT_RADIX;
- ptrdiff_t ddig = (dbits + LOG2_FLT_RADIX - 1) / LOG2_FLT_RADIX;
+ ptrdiff_t ndig = mpz_sizeinbase (*n, FLT_RADIX);
+ ptrdiff_t ddig = mpz_sizeinbase (*d, FLT_RADIX);
/* Scale with SCALE when doing integer division. That is, compute
(N * FLT_RADIX**SCALE) / D [or, if SCALE is negative, N / (D *
FLT_RADIX**-SCALE)] as a bignum, convert the bignum to double,
then divide the double by FLT_RADIX**SCALE. */
- ptrdiff_t scale = ddig - ndig + DBL_MANT_DIG + 1;
+ ptrdiff_t scale = ddig - ndig + DBL_MANT_DIG;
if (scale < 0)
{
mpz_mul_2exp (mpz[1], *d, - (scale * LOG2_FLT_RADIX));
@@ -615,7 +611,7 @@ frac_to_double (Lisp_Object numerator, Lisp_Object
denominator)
round to the nearest integer; otherwise, it is less than
FLT_RADIX ** (DBL_MANT_DIG + 1) and round it to the nearest
multiple of FLT_RADIX. Break ties to even. */
- if (mpz_sizeinbase (*q, 2) < DBL_MANT_DIG * LOG2_FLT_RADIX)
+ if (mpz_sizeinbase (*q, FLT_RADIX) <= DBL_MANT_DIG)
{
/* Converting to double will use the whole quotient so add 1 to
its absolute value as per round-to-even; i.e., if the doubled
@@ -746,46 +742,41 @@ decode_time_components (enum timeform form,
ps = ps % 1000000 + 1000000 * (ps % 1000000 < 0);
us = us % 1000000 + 1000000 * (us % 1000000 < 0);
- if (result)
+ Lisp_Object hz;
+ switch (form)
{
- switch (form)
- {
- case TIMEFORM_HI_LO:
- /* Floats and nil were handled above, so it was an integer. */
- mpz_swap (mpz[0], *s);
- result->hz = make_fixnum (1);
- break;
-
- case TIMEFORM_HI_LO_US:
- mpz_set_ui (mpz[0], us);
- mpz_addmul_ui (mpz[0], *s, 1000000);
- result->hz = make_fixnum (1000000);
- break;
-
- case TIMEFORM_HI_LO_US_PS:
- {
- #if FASTER_TIMEFNS && TRILLION <= ULONG_MAX
- unsigned long i = us;
- mpz_set_ui (mpz[0], i * 1000000 + ps);
- mpz_addmul_ui (mpz[0], *s, TRILLION);
- #else
- intmax_t i = us;
- mpz_set_intmax (mpz[0], i * 1000000 + ps);
- mpz_addmul (mpz[0], *s, ztrillion);
- #endif
- result->hz = trillion;
- }
- break;
+ case TIMEFORM_HI_LO:
+ /* Floats and nil were handled above, so it was an integer. */
+ mpz_swap (mpz[0], *s);
+ hz = make_fixnum (1);
+ break;
- default:
- eassume (false);
- }
- result->ticks = make_integer_mpz ();
+ case TIMEFORM_HI_LO_US:
+ mpz_set_ui (mpz[0], us);
+ mpz_addmul_ui (mpz[0], *s, 1000000);
+ hz = make_fixnum (1000000);
+ break;
+
+ case TIMEFORM_HI_LO_US_PS:
+ {
+ #if FASTER_TIMEFNS && TRILLION <= ULONG_MAX
+ unsigned long i = us;
+ mpz_set_ui (mpz[0], i * 1000000 + ps);
+ mpz_addmul_ui (mpz[0], *s, TRILLION);
+ #else
+ intmax_t i = us;
+ mpz_set_intmax (mpz[0], i * 1000000 + ps);
+ mpz_addmul (mpz[0], *s, ztrillion);
+ #endif
+ hz = trillion;
+ }
+ break;
+
+ default:
+ eassume (false);
}
- else
- *dresult = mpz_get_d (*s) + (us * 1e6L + ps) / 1e12L;
- return 0;
+ return decode_ticks_hz (make_integer_mpz (), hz, result, dresult);
}
enum { DECODE_SECS_ONLY = WARN_OBSOLETE_TIMESTAMPS + 1 };
diff --git a/test/src/timefns-tests.el b/test/src/timefns-tests.el
index 3967590..b24d443 100644
--- a/test/src/timefns-tests.el
+++ b/test/src/timefns-tests.el
@@ -220,6 +220,9 @@ a fixed place on the right and are padded on the left."
'(23752 27217))))
(ert-deftest float-time-precision ()
+ (should (= (float-time '(0 1 0 4025)) 1.000000004025))
+ (should (= (float-time '(1000000004025 . 1000000000000)) 1.000000004025))
+
(should (< 0 (float-time '(1 . 10000000000))))
(should (< (float-time '(-1 . 10000000000)) 0))