libtom / libtommath

LibTomMath is a free open source portable number theoretic multiple-precision integer library written entirely in C.
https://www.libtom.net
Other
650 stars 194 forks source link

Unexpected behavior of `mp_prime_next_prime()` #333

Closed sjaeckel closed 5 years ago

sjaeckel commented 5 years ago

I'll just leave this here... :-D

the diff is based on #332

$ git diff -U1
diff --git a/demo/test.c b/demo/test.c
index 38a75a9..5832472 100644
--- a/demo/test.c
+++ b/demo/test.c
@@ -1135,3 +1135,3 @@ static int test_mp_read_radix(void)

-   while (0) {
+   while (01) {
       char *s = fgets(buf, sizeof(buf), stdin);
@@ -1139,2 +1139,4 @@ static int test_mp_read_radix(void)
       mp_read_radix(&a, buf, 10);
+      mp_to_radix(&a, buf, sizeof(buf), 10);
+      printf("%s, %lu\n", buf, a.dp[0] & 3);
       mp_prime_next_prime(&a, 5, 1);
$ make test_standalone -j9 >/dev/null 
$ ./test read_radix
Digit size 64 Bit 
Size of mp_digit: 8
Size of mp_word: 16
MP_DIGIT_BIT: 60
MP_PREC: 32
TEST mp_read_radix

a == 56
a == 456
a == 123456
8
8, 0
7, 3
czurnieden commented 5 years ago

And I thought I could relax for a day or two and enjoy the last days of summer? *sigh* Nah, you're right, that would be boring! ;-)

+ while (01) {

Did it come with greetings from Rear Admiral G. Hopper? ;-)

Could it be a problem with the specific diet you fed char *s = fgets(buf, sizeof(buf), stdin); with?

But mp_prime_next_prime should return an error in this case, will look into it.

sjaeckel commented 5 years ago

And I thought I could relax for a day or two and enjoy the last days of summer?

please do so!

Did it come with greetings from Rear Admiral G. Hopper? ;-)

of course it did! sir! ;-)

Could it be a problem with the specific diet you fed char *s = fgets(buf, sizeof(buf), stdin); with?

But mp_prime_next_prime should return an error in this case, will look into it.

uhm, why should it return an error?

I gave an 8 to mp_prime_next_prime and I wouldn't expect to get 7 as return value as that's obviously the previous prime or did I understand 'next prime' wrong?

sjaeckel commented 5 years ago

Is it possible that this is the fix?

$ git diff -U1 --ignore-space-change 
diff --git a/bn_mp_prime_next_prime.c b/bn_mp_prime_next_prime.c
index aaa821b..94e2efa 100644
--- a/bn_mp_prime_next_prime.c
+++ b/bn_mp_prime_next_prime.c
@@ -33,3 +33,2 @@ mp_err mp_prime_next_prime(mp_int *a, int t, int bbs_style)
                 */
-               if ((s_mp_prime_tab[x + 1] & 3u) != 3u) {
                /* scan upwards for a prime congruent to 3 mod 4 */
@@ -41,3 +40,2 @@ mp_err mp_prime_next_prime(mp_int *a, int t, int bbs_style)
                }
-               }
             } else {
czurnieden commented 5 years ago

uhm, why should it return an error?

Was just a wild guess without even taking a quick glimpe at the code of bn_mp_prime_next_prime.c. My bad, sorry.

did I understand 'next prime' wrong?

Such a small number gets handled by the table…uh, you have found it.

Is it possible that this is the fix?

No, it shouldn't. If the outer loop works correctly, we should enter that check with x pointing to 7 and x+1 pointing to 11 which is congruent to 3 mod 4, so the output should be 11. But there is no else to that if, that is, if x+1 is actually the next prime congruent to 3 mod 4 it does not return 11 but goes into the inner loop. The next one would be 19 and that should be returned if I am not completely mistaken. So why 7? *scratches head*

No, wait. Shouldn't if (mp_cmp_d(a, s_mp_prime_tab[x]) != MP_LT) { read if (mp_cmp_d(a, s_mp_prime_tab[x]) != MP_GT) { ?

I mean, we look for a number smaller than or equal which means it must not be greater than, that is in LTM speak: != MP_GT Or am I off the rail again?

But that's how far I can get with just staring at it. If that isn't it I need to open the large toolbox and step through it with a debugger which might take a moment longer because my machine is in parts for the seminally de-dusting (and a cap in the PSU looks a bit bulgy and needs replacement, too).

czurnieden commented 5 years ago

Yes, seems to have been the lack of the else in the loop. Please check https://github.com/libtom/libtommath/pull/334

czurnieden commented 5 years ago

There is no standalone test for mp_prime_next_prime() itself. Add one?