diff --git a/bn_mp_prime_strong_lucas_selfridge.c b/bn_mp_prime_strong_lucas_selfridge.c index ca15e3c..5a94f8e 100644 --- a/bn_mp_prime_strong_lucas_selfridge.c +++ b/bn_mp_prime_strong_lucas_selfridge.c @@ -58,9 +58,6 @@ LBL_MPMULSI_ERR: mp_clear(&t); return err; } - - - /* Strong Lucas-Selfridge test. returns MP_YES if it is a strong L-S prime, MP_NO if it is composite @@ -85,10 +82,9 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) /* CZ TODO: Some of them need the full 32 bit, hence the (temporary) exclusion of MP_8BIT */ int32_t D, Ds, J, sign, P, Q, r, s, u, Nbits; int e; - int isset; + int isset, oddness; *result = MP_NO; - /* Find the first element D in the sequence {5, -7, 9, -11, 13, ...} such that Jacobi(D,N) = -1 (Selfridge's algorithm). Theory @@ -138,6 +134,8 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) } } + + P = 1; /* Selfridge's choice */ Q = (1 - Ds) / 4; /* Required so D = P*P - 4*Q */ @@ -232,6 +230,7 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) } Nbits = mp_count_bits(&Dz); + for (u = 1; u < Nbits; u++) { /* zero bit off, already accounted for */ /* Formulas for doubling of indices (carried out mod N). Note that * the indices denoted as "2m" are actually powers of 2, specifically @@ -267,7 +266,6 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) if ((e = mp_mul_2(&Qmz, &Q2mz)) != MP_OKAY) { goto LBL_LS_ERR; } - if ((isset = mp_get_bit(&Dz, u)) == MP_VAL) { e = isset; goto LBL_LS_ERR; @@ -280,7 +278,6 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) * * Be careful with division by 2 (mod N)! */ - if ((e = mp_mul(&U2mz, &Vz, &T1z)) != MP_OKAY) { goto LBL_LS_ERR; } @@ -309,10 +306,11 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) * Thomas R. Nicely used GMP's mpz_fdiv_q_2exp(). * But mp_div_2() does not do so, it is truncating instead. */ + oddness = mp_isodd(&Uz); if ((e = mp_div_2(&Uz, &Uz)) != MP_OKAY) { goto LBL_LS_ERR; } - if ((Uz.sign == MP_NEG) && (mp_isodd(&Uz) != MP_NO)) { + if ((Uz.sign == MP_NEG) && (oddness != MP_NO)) { if ((e = mp_sub_d(&Uz, 1uL, &Uz)) != MP_OKAY) { goto LBL_LS_ERR; } @@ -325,10 +323,11 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) goto LBL_LS_ERR; } } + oddness = mp_isodd(&Vz); if ((e = mp_div_2(&Vz, &Vz)) != MP_OKAY) { goto LBL_LS_ERR; } - if ((Vz.sign == MP_NEG) && (mp_isodd(&Vz) != MP_NO)) { + if ((Vz.sign == MP_NEG) && (oddness != MP_NO)) { if ((e = mp_sub_d(&Vz, 1uL, &Vz)) != MP_OKAY) { goto LBL_LS_ERR; } diff --git a/demo/demo.c b/demo/demo.c index a93ce1d..0aeae25 100644 --- a/demo/demo.c +++ b/demo/demo.c @@ -715,6 +715,16 @@ int main(void) return EXIT_FAILURE; } } + /* Check regarding problem #143 */ +#ifndef MP_8BIT + mp_read_radix(&a, "FFFFFFFFFFFFFFFFC90FDAA22168C234C4C6628B80DC1CD129024E088A67CC74020BBEA63B139B22514A08798E3404DDEF9519B3CD3A431B302B0A6DF25F14374FE1356D6D51C245E485B576625E7EC6F44C42E9A63A3620FFFFFFFFFFFFFFFF", + 16); + mp_prime_strong_lucas_selfridge(&a, &cnt); + if (cnt != MP_YES) { + printf("\n\nissue #143 - mp_prime_strong_lucas_selfridge FAILED!\n"); + return EXIT_FAILURE; + } +#endif printf("\n\n");