mpfr/mpfr-bernoulli-ziv.patch
2020-06-29 14:29:35 -06:00

144 lines
4.9 KiB
Diff
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

diff -Naurd mpfr-4.0.2-a/PATCHES mpfr-4.0.2-b/PATCHES
--- mpfr-4.0.2-a/PATCHES 2020-06-02 11:26:59.031114881 +0000
+++ mpfr-4.0.2-b/PATCHES 2020-06-02 11:26:59.175114040 +0000
@@ -0,0 +1 @@
+bernoulli-ziv
diff -Naurd mpfr-4.0.2-a/VERSION mpfr-4.0.2-b/VERSION
--- mpfr-4.0.2-a/VERSION 2020-04-03 13:54:03.891945830 +0000
+++ mpfr-4.0.2-b/VERSION 2020-06-02 11:26:59.175114040 +0000
@@ -1 +1 @@
-4.0.2-p7
+4.0.2-p8
diff -Naurd mpfr-4.0.2-a/doc/README.dev mpfr-4.0.2-b/doc/README.dev
--- mpfr-4.0.2-a/doc/README.dev 2019-01-07 14:06:05.000000000 +0000
+++ mpfr-4.0.2-b/doc/README.dev 2020-06-02 11:26:59.055114741 +0000
@@ -540,7 +540,9 @@
By default, a fixed seed is used. Only developers
and testers should change the seed.
-+ MPFR_CHECK_LARGEMEM: Define to enable expensive tests.
++ MPFR_CHECK_LARGEMEM: Define to enable tests that take a lot of memory.
+
++ MPFR_CHECK_EXPENSIVE: Define to enable tests that take a lot of time.
+ MPFR_CHECK_LIBC_PRINTF:
Define to enable comparisons with the printf
diff -Naurd mpfr-4.0.2-a/src/bernoulli.c mpfr-4.0.2-b/src/bernoulli.c
--- mpfr-4.0.2-a/src/bernoulli.c 2019-01-07 13:53:20.000000000 +0000
+++ mpfr-4.0.2-b/src/bernoulli.c 2020-06-02 11:26:59.055114741 +0000
@@ -39,7 +39,7 @@
using Von StaudtClausen theorem, which says that the denominator of B[n]
divides the product of all primes p such that p-1 divides n.
Since B[n] = zeta(n) * 2*n!/(2pi)^n, we compute an approximation of
- d * zeta(n) * 2*n!/(2pi)^n and round it to the nearest integer. */
+ (2n+1)! * zeta(n) * 2*n!/(2pi)^n and round it to the nearest integer. */
static void
mpfr_bernoulli_internal (mpz_t *b, unsigned long n)
{
@@ -86,8 +86,9 @@
mpfr_mul_ui (z, z, n, MPFR_RNDU);
p = mpfr_get_ui (z, MPFR_RNDU); /* (n/e/2/pi)^n <= 2^p */
mpfr_clear (z);
- /* the +14 term ensures no rounding failure up to n=10000 */
- prec += p + mpz_sizeinbase (den, 2) + 14;
+ prec += p + mpz_sizeinbase (den, 2);
+ /* the +2 term ensures no rounding failure up to n=10000 */
+ prec += __gmpfr_ceil_log2 (prec) + 2;
}
try_again:
@@ -184,7 +185,6 @@
mpz_mul_ui (t, t, n + 1);
mpz_divexact (t, t, den); /* t was still n! */
mpz_mul (num, num, t);
- mpz_set_ui (den, 1);
mpfr_clear (y);
mpfr_clear (z);
diff -Naurd mpfr-4.0.2-a/src/mpfr.h mpfr-4.0.2-b/src/mpfr.h
--- mpfr-4.0.2-a/src/mpfr.h 2020-04-03 13:54:03.891945830 +0000
+++ mpfr-4.0.2-b/src/mpfr.h 2020-06-02 11:26:59.175114040 +0000
@@ -27,7 +27,7 @@
#define MPFR_VERSION_MAJOR 4
#define MPFR_VERSION_MINOR 0
#define MPFR_VERSION_PATCHLEVEL 2
-#define MPFR_VERSION_STRING "4.0.2-p7"
+#define MPFR_VERSION_STRING "4.0.2-p8"
/* User macros:
MPFR_USE_FILE: Define it to make MPFR define functions dealing
diff -Naurd mpfr-4.0.2-a/src/version.c mpfr-4.0.2-b/src/version.c
--- mpfr-4.0.2-a/src/version.c 2020-04-03 13:54:03.891945830 +0000
+++ mpfr-4.0.2-b/src/version.c 2020-06-02 11:26:59.175114040 +0000
@@ -25,5 +25,5 @@
const char *
mpfr_get_version (void)
{
- return "4.0.2-p7";
+ return "4.0.2-p8";
}
diff -Naurd mpfr-4.0.2-a/tests/tgamma.c mpfr-4.0.2-b/tests/tgamma.c
--- mpfr-4.0.2-a/tests/tgamma.c 2019-01-07 13:53:20.000000000 +0000
+++ mpfr-4.0.2-b/tests/tgamma.c 2020-06-02 11:26:59.055114741 +0000
@@ -1055,6 +1055,48 @@
set_emax (emax);
}
+/* Bug reported by Frithjof Blomquist on May 19, 2020.
+ For the record, this bug was present since r8981
+ (in mpfr_bernoulli_internal, den was wrongly reset to 1 in case
+ of failure in Ziv's loop). The bug only occurred up from r8986
+ where the initial precision was reduced, but was potentially
+ present in any case of failure of Ziv's loop. */
+static void
+bug20200519 (void)
+{
+ mpfr_prec_t prec = 25093;
+ mpfr_t x, y, z, d;
+ double dd;
+ size_t min_memory_limit, old_memory_limit;
+
+ old_memory_limit = tests_memory_limit;
+ min_memory_limit = 24000000;
+ if (tests_memory_limit > 0 && tests_memory_limit < min_memory_limit)
+ tests_memory_limit = min_memory_limit;
+
+ mpfr_init2 (x, prec);
+ mpfr_init2 (y, prec);
+ mpfr_init2 (z, prec + 100);
+ mpfr_init2 (d, 24);
+ mpfr_set_d (x, 2.5, MPFR_RNDN);
+ mpfr_gamma (y, x, MPFR_RNDN);
+ mpfr_gamma (z, x, MPFR_RNDN);
+ mpfr_sub (d, y, z, MPFR_RNDN);
+ mpfr_mul_2si (d, d, prec - mpfr_get_exp (y), MPFR_RNDN);
+ dd = mpfr_get_d (d, MPFR_RNDN);
+ if (dd < -0.5 || 0.5 < dd)
+ {
+ printf ("Error in bug20200519: dd=%f\n", dd);
+ exit (1);
+ }
+ mpfr_clear (x);
+ mpfr_clear (y);
+ mpfr_clear (z);
+ mpfr_clear (d);
+
+ tests_memory_limit = old_memory_limit;
+}
+
int
main (int argc, char *argv[])
{
@@ -1086,6 +1128,11 @@
data_check ("data/gamma", mpfr_gamma, "mpfr_gamma");
+ /* this test takes about one minute */
+ if (getenv ("MPFR_CHECK_EXPENSIVE") != NULL &&
+ getenv ("MPFR_CHECK_LARGEMEM") != NULL)
+ bug20200519 ();
+
tests_end_mpfr ();
return 0;
}