3
\$\begingroup\$

I'm currently writing a Miller-Rabin primality test using C++ and the GNU MP Bignum Library (GMP). As I've only started learning C++ about a day ago, I'm sure this code is riddled with errors and inefficiencies. I'm looking for several things (in order of importance):

  1. Make the code as fast as possible (really at the expense of anything else)
  2. Make the code more readable
  3. Encapsulating all functions in the miller_rabin() function
  4. Make the length of the code shorter

I'm also looking for any insight into my handling of the GMP code. I know there's a way to reallocate memory but I'm not sure how to do that yet. Any thoughts are greatly appreciated!

Here's the code:

#include <iostream>
#include <gmp.h>
#include <gmpxx.h>

using namespace std;


mpz_class exp_one(mpz_class a, mpz_class b, mpz_class c){
    mpz_class r;
    mpz_powm(r.get_mpz_t(), a.get_mpz_t(), b.get_mpz_t(), c.get_mpz_t());
    return r;
}

mpz_class exp_two(mpz_class a, mpz_class b){
    mpz_class r, c;
    c = 2;
    mpz_powm(r.get_mpz_t(), a.get_mpz_t(), c.get_mpz_t(), b.get_mpz_t());
    return r;
}

mpz_class mr_check_one(mpz_class n){
    if(n % 2 != 0){
        return 0;
    }

    n = n - 1;
    while(n % 2 == 0){
        n = n / 2;
    }
    return n;
}

mpz_class mr_check_two(mpz_class m, mpz_class n){
    mpz_class b, u, j, d;
    u = n - 1;
    j = 0;
    d = 2;

    b = exp_one(d, m, n);

    if(b == u or b == 1){
        return 0;
    }

    while(b != u and j < 7){
        b = exp_two(b, n);
        j++;
    }

    if(b == u){
        return 0;
    } else {
        return 1;
    }
}

mpz_class miller_rabin(mpz_class n){
    mpz_class r;
    r = mr_check_one(n);
    r = mr_check_two(r, n);
    return r;
}


int main(){
    mpz_class n, r;
    n = 2305843009213693951;  //some prime
    r = miller_rabin(n);

    if(r == 0){
        cout << "This number is PRIME!" << endl;
    } else {
        cout << "This number is NOT prime." << endl;
    }
}
\$\endgroup\$
  • 1
    \$\begingroup\$ If you are coding this because you need it, GMP already includes an optimized Miller-Rabin test as mpz_probab_prime_p (gmplib.org/manual/…). If you are doing this as a way of learning C++, ignore this comment :). \$\endgroup\$ – Darhuuk Mar 8 '18 at 16:51
  • \$\begingroup\$ @Darhuuk I honestly had no idea. I was using those docs to figure out the GMP functions and syntax and didn't even come across the number-theoretic functions. After running some tests though, I found that my implementation is actually about 15x faster. I do appreciate that you pointed me in that direction though! \$\endgroup\$ – D. Senack Mar 8 '18 at 22:39
  • \$\begingroup\$ Note that your implementation only does a single iteration of Miller-Rabin for a fixed witness value. I.e. I would expect your code to give a lot of false positives. Have you compared the output of your function v.s. e.g. GMP's one? \$\endgroup\$ – Darhuuk Mar 8 '18 at 22:52
  • \$\begingroup\$ @Darhuuk After some testing, it appears to give no false positives (other than the assumed strong pseudoprimes); it does, however, fail to identify around 3% of primes in a test of all numbers from 2 - 10 billion. The difference in speed is what is truly important to me and I only care that the number is "most likely" prime. I'm really searching for probable primes that I can run against another algorithm. GMP out-performs my implementation with smaller numbers but my implementation tends to take half the time when both are run using numbers with over 10,000 digits. Thanks for the help! \$\endgroup\$ – D. Senack Mar 9 '18 at 3:18
  • \$\begingroup\$ Ok, I see why you don't mind the false positives then. However, Miller-Rabin can not give false negatives (i.e. say a number is not prime when it is), so if you're code is doing that in 3% of the cases, then there's definitely a bug in there. \$\endgroup\$ – Darhuuk Mar 9 '18 at 8:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.