From ee15ef7ec67a8cb3a01487420145045bccbbc772 Mon Sep 17 00:00:00 2001 From: czurnieden Date: Sun, 5 Mar 2023 08:16:59 +0100 Subject: [PATCH] Overhaul of the prime-tests - Removal of the Fermat test mp_prime_fermat - Replacement of the Strong Lucas-Selfridge test with the Extra Strong Lucas test with Robert Baillie's parameters P = 3 and Q = 1 - Additional tests to check the implementations of the Miller-Rabin and Extra Strong Lucas tests --- demo/test.c | 165 +++++++++++++++++++++++++++++++++- doc/bn.tex | 23 ++--- libtommath_VS2008.vcproj | 2 +- makefile | 2 +- makefile.mingw | 2 +- makefile.msvc | 2 +- makefile.shared | 2 +- makefile.unix | 2 +- mp_prime_extra_strong_lucas.c | 140 +++++++++++++++++++++++++++++ mp_prime_fermat.c | 41 --------- mp_prime_is_prime.c | 2 +- sources.cmake | 2 +- tommath.def | 2 +- tommath.h | 10 +-- tommath_class.h | 27 ++++-- 15 files changed, 341 insertions(+), 83 deletions(-) create mode 100644 mp_prime_extra_strong_lucas.c delete mode 100644 mp_prime_fermat.c diff --git a/demo/test.c b/demo/test.c index 16fef5570..662fa7e05 100644 --- a/demo/test.c +++ b/demo/test.c @@ -824,6 +824,165 @@ static int test_mp_prime_rand(void) return EXIT_FAILURE; } +/* Some small pseudoprimes to test the individual implementations */ + +/* Miller-Rabin base 2 */ +static const uint32_t SPSP_2[] = { + 2047, 3277, 4033, 4681, 8321, 15841, 29341, 42799, + 49141, 52633, 65281, 74665, 80581, 85489, 88357, 90751 +}; + +/* Miller-Rabin base 3 */ +static const uint32_t SPSP_3[] = { + 121, 703, 1891, 3281, 8401, 8911, 10585, 12403, 16531, + 18721, 19345, 23521, 31621, 44287, 47197, 55969, 63139, + 74593, 79003, 82513, 87913, 88573, 97567 +}; + +/* SPSP to all bases < 100 */ +/* still needs computing +static const char *SPSP_2_100_LARGE[] = { + "", + "", + "", + "", + "" +}; +*/ + +/* Extra strong Lucas test with Baillie's parameters Q = 1, P = 3 */ +static const uint32_t ESLPSP[] = { + 989, 3239, 5777, 10877, 27971, 29681, 30739, 31631, 39059, 72389, + 73919, 75077, 100127, 113573, 125249, 137549, 137801, 153931, 155819, + 161027, 162133, 189419, 218321, 231703, 249331, 370229, 429479, 430127, + 459191, 473891, 480689, 600059, 621781, 632249, 635627 +}; + +/* + Almost extra strong Lucas test with Baillie's parameters Q = 1, P = 3 + Only those that are not in ESLPSP. + */ +static const uint32_t AESLPSP[] = { + 10469, 154697, 233659, 472453, 629693, 852389, 1091093, 1560437, + 1620673, 1813601, 1969109, 2415739, 2595329, 2756837, 3721549, + 4269341, 5192309, 7045433, 7226669, 7265561 +}; + +/* Some randomly choosen 200 decimal digits large primes (https://primes.utm.edu/lists/small/small2.html) */ +static const char *medium_primes[10] = { + "C8Ckh0vviS3HUPdB1NSrSm+gOodw/f1aQ5+aaH1W6RMB0jVkO6lTaL54O3o7U5BSGUFGxm5gAvisbJamasuLZS8g3ZsJ2JM4Vtn9cQZRfkP6b8V", + "64xDN9FqLBiovZ/9q/EPm0DONpIfn5MbJKHa+IjT0fjAzkg34FpAmad+CwhcpKaiTbZEpErut+DhpVyiQfqBFrgcGnGhhIrMF/XkyY3aVx6E96B", + "8cyuMlENm0vh/eWwgHUpDKqmLyCSsRQZRWvbHpA2jHDZv1EhHkVhceg3OFRZn/aXRBnbdtsc2xO6sWh9KZ5Mo7u9rJgBJMVtDnu094MCExj1YvB", + "BRFZFsYjSz45un8qptnuSqEsy9wV0BzbMpVAB1TrwImENOVIc1cASZNQ/mXG2xtazqgn/juVzFo91XLx9PtIlkcK0L2T6fBNgy8Lc7dSVoKQ+XP", + "Ez/mDl+to2gm69+VdIHI9Q7vaO3DuIdLVT69myM3HYwVBE+G24KffAOUAp3FGrSOU+LtERMiIYIEtxPI7n/DRJtmL2i0+REwGpTMge2d2EpabfB", + "5+Uz1gPFjZJ/nNdEOmOaMouJSGzygo42qz7xOwXn/moSUvBpPjo4twRGbK0+qaeU/RI8yYYxXr3OBP4w+/jgL3mN9GiENDM5LtEKMiQrZ9jIVEb", + "AQ5nD1+G1grv41s/XlK+0YTGyZgr/88PzdQJ8QT9tavisTgyG6k8/80A4HQhnFndskHNAaB2EW5fE7KH3kk7m89s8JnVqkJyGZWSfs1+JlmHLPf", + "3F19vPmM0Ih89KZ04Xmd62QB9F6E2sztT10A7Kcqc44eKvsNHh+JY6Z6gJXkbWg1Iw7xr29QAhEF/o1YAgfutQtpdzHkex06Yd71kPsaZdKXiC5", + "2fIcJ1t/VYCColXGs+ji/txNMEXn2FXdowLzlo7QKqzAWHdAbwtltSO5qpSp3OUiEOGUUi3hbyw3iQRE8nFJaikJ89Wdox6vpPtIsc3QRjexMnv", + "8aOicQ5gIbFCarFUgSgzh40LpuZ0jjK1u48/YT+C0h1dAQ8CIEgZjHZT+5/7cCRGmJlo+XCp7S41MSQ2ZNRSJh2texRYtvAXBAZfR8A8twl316P" +}; + +#define ARR_LENGTH(a) ((int)(sizeof((a))/sizeof((a)[0]))) + +static int test_mp_prime_miller_rabin(void) +{ + mp_int a, b, c; + bool result; + int i; + DOR(mp_init_multi(&a, &b, &c, NULL)); + + /* SPSP to base 2 */ + mp_set(&b, 2u); + for (i = 0; i < ARR_LENGTH(SPSP_2); i++) { + result = false; + mp_set_u32(&a, SPSP_2[i]); + DO(mp_prime_miller_rabin(&a, &b, &result)); + EXPECT(result == true); + } + + /* Some larger primes to check for false negatives */ + for (i = 0; i < 10; i++) { + result = false; + DO(mp_read_radix(&a, medium_primes[i], 64)); + DO(mp_prime_miller_rabin(&a, &b, &result)); + EXPECT(result == true); + } + /* Some semi-primes */ + for (i = 0; i < 5; i += 2) { + result = false; + DO(mp_read_radix(&a, medium_primes[i], 64)); + DO(mp_read_radix(&c, medium_primes[i+1], 64)); + DO(mp_mul(&a, &c, &a)); + DO(mp_prime_miller_rabin(&a, &b, &result)); + EXPECT(result == false); + } + + /* SPSP to base 3 */ + mp_set(&b, 3u); + for (i = 0; i < ARR_LENGTH(SPSP_3); i++) { + result = false; + mp_set_u32(&a, SPSP_3[i]); + DO(mp_prime_miller_rabin(&a, &b, &result)); + EXPECT(result == true); + } + + mp_clear_multi(&a, &b, &c, NULL); + return EXIT_SUCCESS; +LBL_ERR: + mp_clear_multi(&a, &b, &c, NULL); + return EXIT_FAILURE; +} + + +static int test_mp_prime_extra_strong_lucas(void) +{ + mp_int a, b; + bool result; + int i; + + DOR(mp_init_multi(&a, &b, NULL)); + + /* Check Extra Strong pseudoprimes */ + for (i = 0; i < ARR_LENGTH(ESLPSP); i++) { + result = false; + mp_set_u32(&a, ESLPSP[i]); + DO(mp_prime_extra_strong_lucas(&a, &result)); + EXPECT(result == true); + } + + /* Check Almost Extra Strong pseudoprimes (not in ESLPSP) */ + for (i = 0; i < ARR_LENGTH(AESLPSP); i++) { + result = false; + mp_set_u32(&a, AESLPSP[i]); + DO(mp_prime_extra_strong_lucas(&a, &result)); + EXPECT(result == false); + } + + /* Some larger primes to check for false negatives */ + for (i = 0; i < 10; i++) { + result = false; + DO(mp_read_radix(&a, medium_primes[i], 64)); + DO(mp_prime_extra_strong_lucas(&a, &result)); + EXPECT(result == true); + } + + /* Some semi-primes */ + for (i = 0; i < 5; i++) { + result = false; + DO(mp_read_radix(&a, medium_primes[i], 64)); + DO(mp_read_radix(&a, medium_primes[i+1], 64)); + DO(mp_mul(&a, &b, &a)); + DO(mp_prime_extra_strong_lucas(&a, &result)); + EXPECT(result == false); + } + + mp_clear_multi(&a, &b, NULL); + return EXIT_SUCCESS; +LBL_ERR: + mp_clear_multi(&a, &b, NULL); + return EXIT_FAILURE; +} + static int test_mp_prime_is_prime(void) { int ix; @@ -893,7 +1052,7 @@ static int test_mp_prime_is_prime(void) DO(mp_read_radix(&a, "FFFFFFFFFFFFFFFFC90FDAA22168C234C4C6628B80DC1CD129024E088A67CC74020BBEA63B139B22514A08798E3404DDEF9519B3CD3A431B302B0A6DF25F14374FE1356D6D51C245E485B576625E7EC6F44C42E9A63A3620FFFFFFFFFFFFFFFF", 16)); - DO(mp_prime_strong_lucas_selfridge(&a, &cnt)); + DO(mp_prime_extra_strong_lucas(&a, &cnt)); /* large problem */ EXPECT(cnt); if ((e != MP_OKAY) || !cnt) { @@ -1358,7 +1517,7 @@ static int test_mp_log_n(void) mp_int a; mp_digit d; int base, lb, size; - const int max_base = MP_MIN(INT_MAX, MP_DIGIT_MAX); + const mp_digit max_base = MP_MIN(INT_MAX, MP_DIGIT_MAX); DOR(mp_init(&a)); @@ -2236,6 +2395,8 @@ static int unit_tests(int argc, char **argv) T1(mp_montgomery_reduce, MP_MONTGOMERY_REDUCE), T1(mp_root_n, MP_ROOT_N), T1(mp_or, MP_OR), + T1(mp_prime_extra_strong_lucas, MP_PRIME_EXTRA_STRONG_LUCAS), + T1(mp_prime_miller_rabin, MP_PRIME_MILLER_RABIN), T1(mp_prime_is_prime, MP_PRIME_IS_PRIME), T1(mp_prime_next_prime, MP_PRIME_NEXT_PRIME), T1(mp_prime_rand, MP_PRIME_RAND), diff --git a/doc/bn.tex b/doc/bn.tex index 566b3be32..5cc7bc62c 100644 --- a/doc/bn.tex +++ b/doc/bn.tex @@ -2014,15 +2014,6 @@ \subsection{Example} \chapter{Prime Numbers} -\section{Fermat Test} -\index{mp\_prime\_fermat} -\begin{alltt} -mp_err mp_prime_fermat (const mp_int *a, const mp_int *b, int *result) -\end{alltt} -Performs a Fermat primality test to the base $b$. That is it computes $b^a \mbox{ mod }a$ and -tests whether the value is equal to $b$ or not. If the values are equal then $a$ is probably prime -and $result$ is set to one. Otherwise $result$ is set to zero. - \section{Miller--Rabin Test} \index{mp\_prime\_miller\_rabin} \begin{alltt} @@ -2032,9 +2023,6 @@ \section{Miller--Rabin Test} test and is very hard to fool (besides with Carmichael numbers). If $a$ passes the test (therefore is probably prime) $result$ is set to one. Otherwise $result$ is set to zero. -Note that it is suggested that you use the Miller--Rabin test instead of the Fermat test since all -of the failures of Miller--Rabin are a subset of the failures of the Fermat test. - \subsection{Required Number of Tests} Generally to ensure a number is very likely to be prime you have to perform the Miller--Rabin with at least a half--dozen or so unique bases. However, it has been proven that the probability of @@ -2229,12 +2217,12 @@ \subsection{Required Number of Tests} See also table C.1 in FIPS 186-4. -\section{Strong Lucas--Selfridge Test} -\index{mp\_prime\_strong\_lucas\_selfridge} +\section{Extra Strong Lucas Test} +\index{mp\_prime\_extra\_strong\_lucas} \begin{alltt} -mp_err mp_prime_strong_lucas_selfridge(const mp_int *a, bool *result) +mp_err mp_prime_extra_strong_lucas(const mp_int *a, bool *result) \end{alltt} -Performs a strong Lucas--Selfridge test. The strong Lucas--Selfridge test together with the +Performs a extra strong Lucas test. The extra strong Lucas test together with the Rabin--Miller test with bases $2$ and $3$ resemble the BPSW test. The single internal use is a compile--time option in \texttt{mp\_prime\_is\_prime} and can be excluded from the Libtommath build if not needed. @@ -2246,8 +2234,7 @@ \section{Frobenius (Underwood) Test} \end{alltt} Performs the variant of the Frobenius test as described by Paul Underwood. It can be included at build--time if the preprocessor macro \texttt{LTM\_USE\_FROBENIUS\_TEST} is defined and will be -used -instead of the Lucas--Selfridge test. +used instead of the extra strong Lucas test. It returns \texttt{MP\_ITER} if the number of iterations is exhausted, assumes a composite as the input and sets \texttt{result} accordingly. This will reduce the set of available pseudoprimes by a diff --git a/libtommath_VS2008.vcproj b/libtommath_VS2008.vcproj index e27aa9860..a39776cd9 100644 --- a/libtommath_VS2008.vcproj +++ b/libtommath_VS2008.vcproj @@ -609,7 +609,7 @@ > = 0; i--) { + if ((err = mp_mul(&Vk, &Vk1, &tmp1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_add(&tmp1, N, &tmp1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_sub_d(&tmp1, P, &tmp1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_mod(&tmp1, N, &tmp1)) != MP_OKAY) goto LTM_ERR; + if (s_mp_get_bit(&s, i) == 1) { + if ((err = mp_copy(&tmp1, &Vk)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_sqr(&Vk1, &Vk1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_add(&Vk1, &Nm2, &Vk1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_mod(&Vk1, N, &Vk1)) != MP_OKAY) goto LTM_ERR; + } else { + if ((err = mp_copy(&tmp1, &Vk1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_sqr(&Vk, &Vk)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_add(&Vk, &Nm2, &Vk)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_mod(&Vk, N, &Vk)) != MP_OKAY) goto LTM_ERR; + } + } + /* Check Vk = \pm 2 mod N */ + if ((mp_cmp_d(&Vk, 2u) == MP_EQ) || (mp_cmp(&Vk, &Nm2) == MP_EQ)) { + /* Compute Uk as described in [2] and check for Uk == 0 */ + if ((err = mp_mul_d(&Vk, P, &tmp1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_mul_2(&Vk1, &Vk1)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_sub(&tmp1, &Vk1, &tmp1)) != MP_OKAY) goto LTM_ERR; + tmp1.sign = MP_ZPOS; + if ((err = mp_mod(&tmp1, N, &Vk1)) != MP_OKAY) goto LTM_ERR; + if (mp_iszero(&Vk1)) { + *result = true; + goto LTM_ERR; + } + } + /* Check for V_{2^t k} = 0 mod N for some t 0 <= t < r - 1 ([3])*/ + for (i = 0; i < (r - 1); i++) { + if (mp_iszero(&Vk)) { + *result = true; + goto LTM_ERR; + } + /* Loop: V(2k) = V(k)^2 - 2 and four minus two is two again */ + if (mp_cmp_d(&Vk, 2u) == MP_EQ) { + goto LTM_ERR; + } + if ((err = mp_sqr(&Vk, &Vk)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_sub_d(&Vk, 2u, &Vk)) != MP_OKAY) goto LTM_ERR; + if ((err = mp_mod(&Vk, N, &Vk)) != MP_OKAY) goto LTM_ERR; + } + +LTM_ERR: + mp_clear_multi(&Dz, &s, &Vk, &Vk1, &tmp1, &Nm2, NULL); + return err; +} +#endif + + + +#endif diff --git a/mp_prime_fermat.c b/mp_prime_fermat.c deleted file mode 100644 index ac8116fef..000000000 --- a/mp_prime_fermat.c +++ /dev/null @@ -1,41 +0,0 @@ -#include "tommath_private.h" -#ifdef MP_PRIME_FERMAT_C -/* LibTomMath, multiple-precision integer library -- Tom St Denis */ -/* SPDX-License-Identifier: Unlicense */ - -/* performs one Fermat test. - * - * If "a" were prime then b**a == b (mod a) since the order of - * the multiplicative sub-group would be phi(a) = a-1. That means - * it would be the same as b**(a mod (a-1)) == b**1 == b (mod a). - * - * Sets result to 1 if the congruence holds, or zero otherwise. - */ -mp_err mp_prime_fermat(const mp_int *a, const mp_int *b, bool *result) -{ - mp_int t; - mp_err err; - - /* ensure b > 1 */ - if (mp_cmp_d(b, 1uL) != MP_GT) { - return MP_VAL; - } - - /* init t */ - if ((err = mp_init(&t)) != MP_OKAY) { - return err; - } - - /* compute t = b**a mod a */ - if ((err = mp_exptmod(b, a, a, &t)) != MP_OKAY) { - goto LBL_ERR; - } - - /* is it equal to b? */ - *result = mp_cmp(&t, b) == MP_EQ; - -LBL_ERR: - mp_clear(&t); - return err; -} -#endif diff --git a/mp_prime_is_prime.c b/mp_prime_is_prime.c index bb24f5944..8e78a5497 100644 --- a/mp_prime_is_prime.c +++ b/mp_prime_is_prime.c @@ -107,7 +107,7 @@ mp_err mp_prime_is_prime(const mp_int *a, int t, bool *result) goto LBL_B; } #else - if ((err = mp_prime_strong_lucas_selfridge(a, &res)) != MP_OKAY) { + if ((err = mp_prime_extra_strong_lucas(a, &res)) != MP_OKAY) { goto LBL_B; } if (!res) { diff --git a/sources.cmake b/sources.cmake index 797d461e2..465f56946 100644 --- a/sources.cmake +++ b/sources.cmake @@ -76,7 +76,7 @@ mp_neg.c mp_or.c mp_pack.c mp_pack_count.c -mp_prime_fermat.c +mp_prime_extra_strong_lucas.c mp_prime_frobenius_underwood.c mp_prime_is_prime.c mp_prime_miller_rabin.c diff --git a/tommath.def b/tommath.def index 1af4e9e79..1f79f0621 100644 --- a/tommath.def +++ b/tommath.def @@ -79,7 +79,7 @@ EXPORTS mp_or mp_pack mp_pack_count - mp_prime_fermat + mp_prime_extra_strong_lucas mp_prime_frobenius_underwood mp_prime_is_prime mp_prime_miller_rabin diff --git a/tommath.h b/tommath.h index 44c92b22a..ae2e4470e 100644 --- a/tommath.h +++ b/tommath.h @@ -499,11 +499,6 @@ mp_err mp_hash(mp_int *a, mp_hval *hash) MP_WUR; /* ---> Primes <--- */ -/* performs one Fermat test of "a" using base "b". - * Sets result to 0 if composite or 1 if probable prime - */ -mp_err mp_prime_fermat(const mp_int *a, const mp_int *b, bool *result) MP_WUR; - /* performs one Miller-Rabin test of "a" using base "b". * Sets result to 0 if composite or 1 if probable prime */ @@ -514,10 +509,11 @@ mp_err mp_prime_miller_rabin(const mp_int *a, const mp_int *b, bool *result) MP_ */ int mp_prime_rabin_miller_trials(int size) MP_WUR; -/* performs one strong Lucas-Selfridge test of "a". + +/* performs one extra strong Lucas-Selfridge test of "a". * Sets result to 0 if composite or 1 if probable prime */ -mp_err mp_prime_strong_lucas_selfridge(const mp_int *a, bool *result) MP_WUR; +mp_err mp_prime_extra_strong_lucas(const mp_int *a, bool *result) MP_WUR; /* performs one Frobenius test of "a" as described by Paul Underwood. * Sets result to 0 if composite or 1 if probable prime diff --git a/tommath_class.h b/tommath_class.h index becbcb193..84599b07f 100644 --- a/tommath_class.h +++ b/tommath_class.h @@ -85,7 +85,7 @@ # define MP_OR_C # define MP_PACK_C # define MP_PACK_COUNT_C -# define MP_PRIME_FERMAT_C +# define MP_PRIME_EXTRA_STRONG_LUCAS_C # define MP_PRIME_FROBENIUS_UNDERWOOD_C # define MP_PRIME_IS_PRIME_C # define MP_PRIME_MILLER_RABIN_C @@ -582,12 +582,27 @@ # define MP_COUNT_BITS_C #endif -#if defined(MP_PRIME_FERMAT_C) -# define MP_CLEAR_C +#if defined(MP_PRIME_EXTRA_STRONG_LUCAS_C) +# define MP_ADD_C +# define MP_ADD_D_C +# define MP_CLEAR_MULTI_C # define MP_CMP_C # define MP_CMP_D_C -# define MP_EXPTMOD_C -# define MP_INIT_C +# define MP_CNT_LSB_C +# define MP_COPY_C +# define MP_COUNT_BITS_C +# define MP_DIV_2D_C +# define MP_INIT_MULTI_C +# define MP_KRONECKER_C +# define MP_MOD_C +# define MP_MUL_2_C +# define MP_MUL_C +# define MP_MUL_D_C +# define MP_SET_C +# define MP_SET_U32_C +# define MP_SUB_C +# define MP_SUB_D_C +# define S_MP_GET_BIT_C #endif #if defined(MP_PRIME_FROBENIUS_UNDERWOOD_C) @@ -619,8 +634,8 @@ # define MP_DIV_2D_C # define MP_INIT_SET_C # define MP_IS_SQUARE_C +# define MP_PRIME_EXTRA_STRONG_LUCAS_C # define MP_PRIME_MILLER_RABIN_C -# define MP_PRIME_STRONG_LUCAS_SELFRIDGE_C # define MP_RAND_C # define MP_READ_RADIX_C # define MP_SET_C