# Miller-Rabin prime test

I have an implementation of the Miller-Rabin prime test, coded in Python 3.

I want to use that as a comparison to a prime test which I coded myself. Unfortunately I have no idea how "fast" the Miller-Rabin test should be on my machine.

With my code I can test the number n = 2¹⁰⁰⁰⁰+1 in less than 4 seconds. But I don't know whether this is a good time.

The implemented version of gmpy2 is much much faster. But I think it's not a good comparison because my own test is written in python and gmpy2 is written in C, I suppose.

import gmpy2
import random
import time

def rab_test(N, rounds = 40):
if N in [2, 3]:
return True
if N % 2 == 0:
return False
n = N // 2
t = 1
while n % 2 == 0:
n = n // 2
t += 1
for _ in range(0, rounds):
a = random.randrange(2, N-1)
u = gmpy2.powmod(a, n, N)
if u == 1 or N-u == 1:
return True
for k in range(1, t):
u = gmpy2.powmod(u, 2, N)
if N-u == 1:
return True
return False


I tested 2²⁰⁰⁰+1 and got the result after 0.5 seconds. For 10²⁰⁰⁰+1, it needed 4.8 seconds

start = time.time()
print(rab_test(2**2000+1))
print(time.time() - start)

start = time.time()
print(rab_test(10**2000+1))
print(time.time() - start)


For (2⁵⁰⁰¹+1)*2¹⁰⁰⁰⁰+1, it needed exactly 60 seconds.

• If you want to use it as the baseline to compare another algorithm, then the performance relative to a C implementation isn't important, I think - only how well your own tests scales against it. Mar 12 at 15:37

As written, the OP code gives me little reason to believe it is a faithful implementation of Miller's paper nor of wikipedia's retelling of it.

Using an opaque identifier like u can be very helpful to someone seeking to understand your implementation, when it matches a variable from your cited reference. You cited wikipedia, yet invented new names unrelated to the variable names used there. This fails to inspire confidence in the implementation.

# wrong name

def rab_test might be a good way to introduce the analysis of a rabies test, but in this context it is less helpful than e.g. is_a_miller_rabin_prime(). As stated it's unclear how we should interpret a True result until we examine some special cases in the code.

    n = N // 2


That would be an unfortunate choice of identifiers even if some cited author used two different N's to represent two different quantities. The pep-8 advice on choosing names is sometimes at variance with long established mathematical conventions. Here, I would ask you to please at least choose names that we could recite over the telephone when discussing details with a colleague, without resorting to verbose descriptions of "no not that N, the other N in Fraktur font" or whatever the original naming motivation was.

# docstring

The OP source code has no comments of any kind. A """docstring""" would be an excellent place to offer one or more literature citations.

# use tqdm

Performance is of interest for this code, which can take some seconds to produce a result. You may find it convenient to instrument the outer loop with a progress bar:

from tqdm import tqdm
...
for _ in tqdm(range(0, rounds)):


As you try different input values and code variations this lets you quickly assess the cost of each round.

Consider using numpy to compute the gmpy2.powmod(u, 2, N) recurrence relation as a single vector of t objects, in the interest of speed.

## Documentation

You could add documentation at the top of your file to state the purpose of the code, for example:

'''
The Miller-Rabin prime test is an algorithm which determines
whether a given number is likely to be prime (from Wikipedia).
'''


Your function could have a docstring which describes valid inputs and the return value.

## Lint check

pylint identified a few style issues.

Variable k is unused in this line:

    for k in range(1, t):


You can use _ as you did in your other for loop.

## Naming

It is common for functions that return a boolean value to have an is_ prefix. For example, instead of rab_test, consider naming the function something like is_likely_prime.

Most of the variable names are not very meaningful (a, t, u, etc.). Consider giving them longer, more meaningful names, or at least add comments describing what they refer to. If they are based on an equation, you can add a reference to some external documentation.

## Simple

This could be simplified:

    n = n // 2


as:

    n //= 2


## Performance

As you mentioned, the gmpy2 library your code uses is implemented in C. This means you are already potentially taking advantage of fast code.

## Test code

I added a loop around your calls to make it easier for others to add more examples:

for n in [2**2000+1, 10**2000+1]:
start = time.time()
print(rab_test(n))
print(time.time() - start)

• Thank you, there are helpful hints and I will take note of them. Nevertheless, my actual question is not answered. I was wondering whether it is meaningful to compare my test (written in python) to an implementation of the miller-rabin written in C. Mar 10 at 17:22
• @Lereu: You're welcome. Since you are a new user, another way to say "thanks" for those who answered your question is to upvote/accept Mar 15 at 9:26