Torbjorn Granlund wrote: > The very old factoring code cut from an now obsolete version GMP does > not pass proper arguments to the mpz_probab_prime_p function. It ask > for 3 Miller-Rabin tests only, which is not sufficient.
Hi Torbjorn Thank you for the patch and explanation. I've converted that into the commit below in your name. Please proofread it and let me know if you'd like to change anything. I tweaked the patch to change MR_REPS from a #define to an enum and to add the comment just preceding. I'll add NEWS and tests separately. >From ea6dd126e6452504f9fa1d6708d25473e2c27e67 Mon Sep 17 00:00:00 2001 From: Torbjorn Granlund <t...@gmplib.org> Date: Tue, 4 Sep 2012 16:22:47 +0200 Subject: [PATCH] factor: don't ever declare composites to be prime The multiple-precision factoring code (with HAVE_GMP) was copied from a now-obsolete version of GMP that did not pass proper arguments to the mpz_probab_prime_p function. It makes that code perform no more than 3 Miller-Rabin tests only, which is not sufficient. A Miller-Rabin test will detect composites with at least a probability of 3/4. For a uniform random composite, the probability will actually by much higher. Or put another way, of the N-3 possible Miller-Rabin tests for checking the composite N, there is no number N for which more than (N-3)/4 of the tests will fail to detect the number as a composite. For most numbers N the number of "false witnesses" will be much, much lower. Problem numbers are of the for N=pq, p,q prime and (p-1)/(q-1) = s, where s is a small integer. (There are other problem forms too, incvolving 3 or more prime factors.) When s = 2, we get the 3/4 factor. It is easy to find numbers of that form that cause coreutils' factor to fail: 465658903 2242724851 6635692801 17709149503 17754345703 20889169003 42743470771 54890944111 72047131003 85862644003 98275842811 114654168091 117225546301 ... There are 9008992 composites of the form with s=2 below 2^64. With 3 Miller-Rabin test, one would expect about 9008992/4^64 = 140766 to be invalidly recognised as primes in that range. * src/factor.c (MR_REPS): Define to 25. (factor_using_pollard_rho): Use MR_REPS, not 3. (print_factors_multi): Likewise. --- src/factor.c | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/factor.c b/src/factor.c index 1d55805..e63e0e0 100644 --- a/src/factor.c +++ b/src/factor.c @@ -153,6 +153,9 @@ factor_using_division (mpz_t t, unsigned int limit) mpz_clear (r); } +/* The number of Miller-Rabin tests we require. */ +enum { MR_REPS = 25 }; + static void factor_using_pollard_rho (mpz_t n, int a_int) { @@ -222,7 +225,7 @@ S4: mpz_div (n, n, g); /* divide by g, before g is overwritten */ - if (!mpz_probab_prime_p (g, 3)) + if (!mpz_probab_prime_p (g, MR_REPS)) { do { @@ -242,7 +245,7 @@ S4: mpz_mod (x, x, n); mpz_mod (x1, x1, n); mpz_mod (y, y, n); - if (mpz_probab_prime_p (n, 3)) + if (mpz_probab_prime_p (n, MR_REPS)) { emit_factor (n); break; @@ -411,7 +414,7 @@ print_factors_multi (mpz_t t) if (mpz_cmp_ui (t, 1) != 0) { debug ("[is number prime?] "); - if (mpz_probab_prime_p (t, 3)) + if (mpz_probab_prime_p (t, MR_REPS)) emit_factor (t); else factor_using_pollard_rho (t, 1); -- 1.7.12.176.g3fc0e4c