0

I'm attempting to debug my implementation of the standard (i.e., probabilistic) Miller-Rabin test. This is the code right now, with a rudimentary test execution in the main method:

#include <iostream>
#include <random>

typedef unsigned long long int ulli;
std::random_device rd;
std::mt19937_64 mt(rd());

ulli upow(ulli x, ulli y) {
    if (y == 0)
        return 1;
    else if (y == 1)
        return x;
    else if (y % 2 == 0)
        return upow(x*x, y/2);
    else
        return x * upow(x*x, (y-1)/2);
}

bool miller_rabin(ulli n, ulli d) {
    std::uniform_int_distribution<ulli> dist(2, n-2);

    ulli a = dist(mt);
    ulli x = upow(a, d) % n;
    if (x == 1 || x == n - 1)
        return true;

    while (d != n - 1) {
        x = x*x % n;
        d *= 2;
        if (x == 1)
            return false;
        if (x == n - 1)
            return true;
    }

    return false;
}

bool is_prime(ulli n, ulli k) {
    if (n <= 1)
        return false;
    if (n <= 3)
        return true;
    if (n % 2 == 0)
        return false;

    ulli d = n - 1;
    while (d % 2 == 0)
        d /= 2;

    for (int i = 0; i < k; ++i) {
        if (!miller_rabin(n, d))
            return false;
    }

    return true;
}


int main() {
    ulli tests = 30;

    std::cout << "is_prime(2): " << is_prime(2, tests) << std::endl;
    std::cout << "is_prime(10): " << is_prime(10, tests) << std::endl;
    std::cout << "is_prime(97): " << is_prime(97, tests) << std::endl;
    std::cout << "is_prime(256): " << is_prime(256, tests) << std::endl;
    std::cout << "is_prime(179424691): " << is_prime(179424691, tests) << std::endl;
    std::cout << "is_prime(32416187563): " << is_prime(32416187563, tests) << std::endl;
}

The problem I'm currently having is this: the last two tests are returning false, which (because they are both prime) should be impossible. Even with as high a number of tests as I'm passing, is_prime still outputs false.

At this point, I'm genuinely lost as to what's wrong. Any ideas?

DaBler
  • 2,695
  • 2
  • 26
  • 46
MrM21632
  • 27
  • 1
  • 8
  • Please learn how to use a debugger, and use it to step through your code. – Rakete1111 Oct 13 '17 at 15:41
  • `% n` is missing inside `upow`, so the `ulli` type overflows. – DaBler Oct 13 '17 at 16:06
  • 1
    @DaBler That seems to have been the issue. It works as expected now, although the last test is still returning false (must still be overflowing). Sorry, by the way, I probably should have used a debugger first. – MrM21632 Oct 13 '17 at 16:48
  • 1
    I have fully solved my problem. I also had to pull out the modular multiplication (e.g., `x = x*x % n`) as well, because that too was overflowing. Now it works exactly as expected. – MrM21632 Oct 13 '17 at 19:18

0 Answers0