changes and bigfixes, see pull-request #113 at https://github.com/libtom/libtommath/pull/113 for details

This commit is contained in:
czurnieden 2018-05-21 22:17:48 +02:00 committed by Steffen Jaeckel
parent 38e8f93bdb
commit 934dd31738
5 changed files with 100 additions and 61 deletions

View File

@ -27,6 +27,21 @@ int mp_get_bit(const mp_int *a, int b)
} }
limb = b / DIGIT_BIT; limb = b / DIGIT_BIT;
/*
* Zero is a special value with the member "used" set to zero.
* Needs to be tested before the check for the upper boundary
* otherwise (limb >= a->used) would be true for a = 0
*/
if(mp_iszero(a)) {
return MP_NO;
}
if (limb >= a->used) {
return MP_VAL;
}
bit = (mp_digit)1 << ((mp_digit)b % DIGIT_BIT); bit = (mp_digit)1 << ((mp_digit)b % DIGIT_BIT);
isset = a->dp[limb] & bit; isset = a->dp[limb] & bit;
return (isset != 0) ? MP_YES : MP_NO; return (isset != 0) ? MP_YES : MP_NO;

View File

@ -101,14 +101,16 @@ int mp_kronecker(const mp_int *a, const mp_int *p, int *c)
} }
if (a1.sign == MP_NEG) { if (a1.sign == MP_NEG) {
// compute k = (-1)^((a1)*(p1-1)/4) * k /*
// a1.dp[0] + 1 cannot overflow because the MSB * Compute k = (-1)^((a1)*(p1-1)/4) * k
// of the type mp_digit is not set by definition * a1.dp[0] + 1 cannot overflow because the MSB
* of the type mp_digit is not set by definition
*/
if ((a1.dp[0] + 1) & p1.dp[0] & 2u) { if ((a1.dp[0] + 1) & p1.dp[0] & 2u) {
k = -k; k = -k;
} }
} else { } else {
// compute k = (-1)^((a1-1)*(p1-1)/4) * k /* compute k = (-1)^((a1-1)*(p1-1)/4) * k */
if (a1.dp[0] & p1.dp[0] & 2u) { if (a1.dp[0] & p1.dp[0] & 2u) {
k = -k; k = -k;
} }
@ -128,9 +130,9 @@ int mp_kronecker(const mp_int *a, const mp_int *p, int *c)
LBL_KRON: LBL_KRON:
mp_clear(&r); mp_clear(&r);
LBL_KRON_0:
mp_clear(&a1);
LBL_KRON_1: LBL_KRON_1:
mp_clear(&a1);
LBL_KRON_0:
mp_clear(&p1); mp_clear(&p1);
return e; return e;
} }

View File

@ -16,17 +16,23 @@
#ifdef MP_8BIT #ifdef MP_8BIT
// floor of positive solution of /*
// (2^16)-1 = (a+4)*(2*a+5) * floor of positive solution of
// TODO: that is too small, would have to use a bigint for a instead * (2^16)-1 = (a+4)*(2*a+5)
* TODO: that is too small, would have to use a bigint for a instead
*/
#define LTM_FROBENIUS_UNDERWOOD_A 177 #define LTM_FROBENIUS_UNDERWOOD_A 177
// commented out to allow Travis's tests to run /*
// Don't forget to switch in back on in production or we'll find it at TDWTF.com! * Commented out to allow Travis's tests to run
//#warning "Frobenius test not fully usable with MP_8BIT!" * Don't forget to switch it back on in production or we'll find it at TDWTF.com!
*/
/* #warning "Frobenius test not fully usable with MP_8BIT!" */
#else #else
// floor of positive solution of /*
// (2^31)-1 = (a+4)*(2*a+5) * floor of positive solution of
// TODO: that might be too small * (2^31)-1 = (a+4)*(2*a+5)
* TODO: that might be too small
*/
#define LTM_FROBENIUS_UNDERWOOD_A 32764 #define LTM_FROBENIUS_UNDERWOOD_A 32764
#endif #endif
int mp_prime_frobenius_underwood(const mp_int *N, int *result) int mp_prime_frobenius_underwood(const mp_int *N, int *result)
@ -43,11 +49,11 @@ int mp_prime_frobenius_underwood(const mp_int *N, int *result)
} }
for (a = 0; a < LTM_FROBENIUS_UNDERWOOD_A; a++) { for (a = 0; a < LTM_FROBENIUS_UNDERWOOD_A; a++) {
//TODO: That's ugly! No, really, it is! /* TODO: That's ugly! No, really, it is! */
if (a==2||a==4||a==7||a==8||a==10||a==14||a==18||a==23||a==26||a==28) { if (a==2||a==4||a==7||a==8||a==10||a==14||a==18||a==23||a==26||a==28) {
continue; continue;
} }
// (32764^2 - 4) < 2^31, no bigint for >MP_8BIT needed) /* (32764^2 - 4) < 2^31, no bigint for >MP_8BIT needed) */
if ((e = mp_set_long(&T1z,(unsigned long)a)) != MP_OKAY) { if ((e = mp_set_long(&T1z,(unsigned long)a)) != MP_OKAY) {
goto LBL_FU_ERR; goto LBL_FU_ERR;
} }
@ -69,7 +75,7 @@ int mp_prime_frobenius_underwood(const mp_int *N, int *result)
} }
if (j == 0) { if (j == 0) {
// composite /* composite */
goto LBL_FU_ERR; goto LBL_FU_ERR;
} }
} }
@ -77,7 +83,7 @@ int mp_prime_frobenius_underwood(const mp_int *N, int *result)
e = MP_VAL; e = MP_VAL;
goto LBL_FU_ERR; goto LBL_FU_ERR;
} }
// Composite if N and (a+4)*(2*a+5) are not coprime /* Composite if N and (a+4)*(2*a+5) are not coprime */
if ((e = mp_set_long(&T1z, (unsigned long)((a+4)*(2*a+5)))) != MP_OKAY) { if ((e = mp_set_long(&T1z, (unsigned long)((a+4)*(2*a+5)))) != MP_OKAY) {
goto LBL_FU_ERR; goto LBL_FU_ERR;
} }
@ -101,14 +107,14 @@ int mp_prime_frobenius_underwood(const mp_int *N, int *result)
for (i = length - 2; i >= 0; i--) { for (i = length - 2; i >= 0; i--) {
/* /*
temp = (sz*(a*sz+2*tz))%N; * temp = (sz*(a*sz+2*tz))%N;
tz = ((tz-sz)*(tz+sz))%N; * tz = ((tz-sz)*(tz+sz))%N;
sz = temp; * sz = temp;
*/ */
if ((e = mp_mul_2(&tz,&T2z)) != MP_OKAY) { if ((e = mp_mul_2(&tz,&T2z)) != MP_OKAY) {
goto LBL_FU_ERR; goto LBL_FU_ERR;
} }
// TODO: is this small saving worth the branch? /* a = 0 at about 50% of the cases (non-square and odd input) */
if (a != 0) { if (a != 0) {
if ((e = mp_mul_d(&sz,(mp_digit)a,&T1z)) != MP_OKAY) { if ((e = mp_mul_d(&sz,(mp_digit)a,&T1z)) != MP_OKAY) {
goto LBL_FU_ERR; goto LBL_FU_ERR;
@ -141,9 +147,9 @@ int mp_prime_frobenius_underwood(const mp_int *N, int *result)
} }
if (isset == MP_YES) { if (isset == MP_YES) {
/* /*
temp = (a+2) * sz + tz * temp = (a+2) * sz + tz
tz = 2 * tz - sz * tz = 2 * tz - sz
sz = temp * sz = temp
*/ */
if (a == 0) { if (a == 0) {
if ((e = mp_mul_2(&sz,&T1z)) != MP_OKAY) { if ((e = mp_mul_2(&sz,&T1z)) != MP_OKAY) {

View File

@ -119,16 +119,16 @@ int mp_prime_is_prime(const mp_int *a, int t, int *result)
t = 8; t = 8;
} }
#else #else
// commented out for testing purposes /* commented out for testing purposes */
//#ifdef LTM_USE_STRONG_LUCAS_SELFRIDGE_TEST /* #ifdef LTM_USE_STRONG_LUCAS_SELFRIDGE_TEST */
if ((err = mp_prime_strong_lucas_selfridge(a, &res)) != MP_OKAY) { if ((err = mp_prime_strong_lucas_selfridge(a, &res)) != MP_OKAY) {
goto LBL_B; goto LBL_B;
} }
if (res == MP_NO) { if (res == MP_NO) {
goto LBL_B; goto LBL_B;
} }
//#endif /* #endif */
// commented out for testing purposes /* commented out for testing purposes */
#ifdef LTM_USE_FROBENIUS_UNDERWOOD_TEST #ifdef LTM_USE_FROBENIUS_UNDERWOOD_TEST
if ((err = mp_prime_frobenius_underwood(a, &res)) != MP_OKAY) { if ((err = mp_prime_frobenius_underwood(a, &res)) != MP_OKAY) {
goto LBL_B; goto LBL_B;
@ -223,12 +223,14 @@ int mp_prime_is_prime(const mp_int *a, int t, int *result)
See Fips 186.4 p. 126ff See Fips 186.4 p. 126ff
*/ */
else if (t > 0) { else if (t > 0) {
// The mp_digit's have a defined bit-size but the size of the /*
// array a.dp is a simple 'int' and this library can not assume full * The mp_digit's have a defined bit-size but the size of the
// compliance to the current C-standard (ISO/IEC 9899:2011) because * array a.dp is a simple 'int' and this library can not assume full
// it gets used for small embeded processors, too. Some of those MCUs * compliance to the current C-standard (ISO/IEC 9899:2011) because
// have compilers that one cannot call standard compliant by any means. * it gets used for small embeded processors, too. Some of those MCUs
// Hence the ugly type-fiddling in the following code. * have compilers that one cannot call standard compliant by any means.
* Hence the ugly type-fiddling in the following code.
*/
size_a = mp_count_bits(a); size_a = mp_count_bits(a);
mask = (1u << floor_ilog2(size_a)) - 1u; mask = (1u << floor_ilog2(size_a)) - 1u;
/* /*
@ -266,16 +268,20 @@ int mp_prime_is_prime(const mp_int *a, int t, int *result)
need to be prime. need to be prime.
*/ */
for (ix = 0; ix < t; ix++) { for (ix = 0; ix < t; ix++) {
// mp_rand() guarantees the first digit to be non-zero /* mp_rand() guarantees the first digit to be non-zero */
if ((err = mp_rand(&b, 1)) != MP_OKAY) { if ((err = mp_rand(&b, 1)) != MP_OKAY) {
goto LBL_B; goto LBL_B;
} }
// Reduce digit before casting because mp_digit might be bigger than /*
// an unsigned int and "mask" on the other side is most probably not. * Reduce digit before casting because mp_digit might be bigger than
* an unsigned int and "mask" on the other side is most probably not.
*/
fips_rand = (unsigned int) (b.dp[0] & (mp_digit) mask); fips_rand = (unsigned int) (b.dp[0] & (mp_digit) mask);
#ifdef MP_8BIT #ifdef MP_8BIT
// One 8-bit digit is too small, so concatenate two if the size of /*
// unsigned int allows for it. * One 8-bit digit is too small, so concatenate two if the size of
* unsigned int allows for it.
*/
if( (sizeof(unsigned int) * CHAR_BIT)/2 >= (sizeof(mp_digit) * CHAR_BIT) ) { if( (sizeof(unsigned int) * CHAR_BIT)/2 >= (sizeof(mp_digit) * CHAR_BIT) ) {
if ((err = mp_rand(&b, 1)) != MP_OKAY) { if ((err = mp_rand(&b, 1)) != MP_OKAY) {
goto LBL_B; goto LBL_B;
@ -285,19 +291,21 @@ int mp_prime_is_prime(const mp_int *a, int t, int *result)
fips_rand &= mask; fips_rand &= mask;
} }
#endif #endif
// Ceil, because small numbers have a right to live, too, /* Ceil, because small numbers have a right to live, too, */
len = (int) ( (fips_rand + DIGIT_BIT) / DIGIT_BIT); len = (int) ( (fips_rand + DIGIT_BIT) / DIGIT_BIT);
// Unlikely. /* Unlikely. */
if(len < 0){ if(len < 0){
ix--; ix--;
continue; continue;
} }
// As mentioned above, one 8-bit digit is too small and /*
// although it can only happen in the unlikely case that * As mentioned above, one 8-bit digit is too small and
// an "unsigned int" is smaller than 16 bit a simple test * although it can only happen in the unlikely case that
// is cheap and the correction even cheaper. * an "unsigned int" is smaller than 16 bit a simple test
* is cheap and the correction even cheaper.
*/
#ifdef MP_8BIT #ifdef MP_8BIT
// All "a" < 2^8 have been caught before /* All "a" < 2^8 have been caught before */
if(len == 1){ if(len == 1){
len++; len++;
} }
@ -305,15 +313,17 @@ int mp_prime_is_prime(const mp_int *a, int t, int *result)
if ((err = mp_rand(&b, len)) != MP_OKAY) { if ((err = mp_rand(&b, len)) != MP_OKAY) {
goto LBL_B; goto LBL_B;
} }
// That number might got too big and the witness has to be /*
// smaller than or equal to "a" * That number might got too big and the witness has to be
* smaller than or equal to "a"
*/
len = mp_count_bits(&b); len = mp_count_bits(&b);
if (len > size_a) { if (len > size_a) {
len = len - size_a; len = len - size_a;
mp_div_2d(&b, len, &b, NULL); mp_div_2d(&b, len, &b, NULL);
} }
// Although the chance for b <= 3 is miniscule, try again. /* Although the chance for b <= 3 is miniscule, try again. */
if (mp_cmp_d(&b,3) != MP_GT) { if (mp_cmp_d(&b,3) != MP_GT) {
ix--; ix--;
continue; continue;

View File

@ -28,16 +28,16 @@
liability arising from its use liability arising from its use
The multi-line comments are made by Thomas R. Nicely and are copied verbatim. The multi-line comments are made by Thomas R. Nicely and are copied verbatim.
Single-line comments are by the code-portist. Additional comments marked "CZ" (without the quotes) are by the code-portist.
(If that name sounds familiar, he is the guy who found the fdiv bug in the (If that name sounds familiar, he is the guy who found the fdiv bug in the
Pentium (P5x, I think) Intel processor) Pentium (P5x, I think) Intel processor)
*/ */
int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result) int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result)
{ {
// TODO: choose better variable names! /* CZ TODO: choose better variable names! */
mp_int Dz, gcd, Np1, Uz, Vz, U2mz, V2mz, Qmz, Q2mz, Qkdz, T1z, T2z, T3z, T4z, Q2kdz; mp_int Dz, gcd, Np1, Uz, Vz, U2mz, V2mz, Qmz, Q2mz, Qkdz, T1z, T2z, T3z, T4z, Q2kdz;
// TODO: Some of them need the full 32 bit, hence the (temporary) exclusion of MP_8BIT /* 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; int32_t D, Ds, J, sign, P, Q, r, s, u, Nbits;
int e = MP_OKAY; int e = MP_OKAY;
int isset; int isset;
@ -131,9 +131,13 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result)
} }
s = mp_cnt_lsb(&Np1); s = mp_cnt_lsb(&Np1);
// this should round towards zero because /* CZ
// Thomas R. Nicely used GMP's mpz_tdiv_q_2exp() * This should round towards zero because
// and mp_div_2d() is equivalent * Thomas R. Nicely used GMP's mpz_tdiv_q_2exp()
* and mp_div_2d() is equivalent. Additionally:
* dividing an even number by two does not produce
* any leftovers.
*/
if ((e = mp_div_2d(&Np1, s, &Dz, NULL)) != MP_OKAY) { if ((e = mp_div_2d(&Np1, s, &Dz, NULL)) != MP_OKAY) {
goto LBL_LS_ERR; goto LBL_LS_ERR;
} }
@ -211,7 +215,7 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result)
if ((e = mp_sqr(&Qmz,&Qmz)) != MP_OKAY) { if ((e = mp_sqr(&Qmz,&Qmz)) != MP_OKAY) {
goto LBL_LS_ERR; goto LBL_LS_ERR;
} }
/* prevents overflow */ // still necessary without a fixed prealloc'd mem.? /* prevents overflow */ /* CZ still necessary without a fixed prealloc'd mem.? */
if ((e = mp_mod(&Qmz,a,&Qmz)) != MP_OKAY) { if ((e = mp_mod(&Qmz,a,&Qmz)) != MP_OKAY) {
goto LBL_LS_ERR; goto LBL_LS_ERR;
} }
@ -255,9 +259,11 @@ int mp_prime_strong_lucas_selfridge(const mp_int *a, int *result)
goto LBL_LS_ERR; goto LBL_LS_ERR;
} }
} }
// This should round towards negative infinity because /* CZ
// Thomas R. Nicely used GMP's mpz_fdiv_q_2exp(). * This should round towards negative infinity because
// But mp_div_2() does not do so, it is truncating instead. * Thomas R. Nicely used GMP's mpz_fdiv_q_2exp().
* But mp_div_2() does not do so, it is truncating instead.
*/
if ((e = mp_div_2(&Uz,&Uz)) != MP_OKAY) { if ((e = mp_div_2(&Uz,&Uz)) != MP_OKAY) {
goto LBL_LS_ERR; goto LBL_LS_ERR;
} }