# Miller-Rabin Primality Test-Witnesses and Liars - Implementing in Python

I have been studying the Miller-Rabin Primality Test, and am interested in implementing a code in Python to count witnesses and liars. In the basic code to determine if a number is probably prime or composite, I would like to incorporate 1) and 2) below to better understand witnesses compared to liars for n values:

$$1)$$ $$a$$ to be tested for all values of $$a$$ from $$1, not random as it is currently coded.

$$2)$$ Then for every $$a$$ in $$1)$$ above, a count of how many of those $$a's$$ are witnesses and how many are non-witnesses(liars)

My ultimate goal is to use this, I'm sure with more modifications to the code, to compare to the Theorem: If n is odd, composite, and n>9, then at least 75% of the elements in $$(\mathbb Z/n\mathbb Z)^x$$ are Miller-Rabin witnesses.

The Python code I am using is as follows:

from random import randrange

def probably_prime(n, k):
"""Return True if n passes k rounds of the Miller-Rabin primality
test (and is probably prime). Return False if n is proved to be
composite.

"""
if n < 2: return False

r, m = 0, n - 1
while m % 2 == 0:
r += 1
m //= 2
for _ in range(k):
a = randrange(2, n - 1)
x = pow(a, m, n)
if x == 1 or x == n - 1:
continue
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
return True

• I am not sure I get your questions. For question 1, do you mean you want to add $\alpha$ as an argument in your algorithm? For question 2, do you mean you want to test every single $\alpha$ from 1 to $n-1$, and find which are Miller-Rabin liars and which are not? – Alexandros Feb 15 at 22:24
• For 1) I would like to test every $a$ in the range, not random as it is now. 2) Yes, that is correct - Thank you – user565684 Feb 15 at 23:06
• Why don't you just run over all $a$ with a for loop, like the one you already use for _? Before the loop on $a$ starts, you can define two variables both set to 0, which you increment later on as needed to count the witnesses & liars. – J.G. Feb 15 at 23:18
• Thank you, I will work on this, but I am very new to coding and do not have a great understanding of it at this point...but I am really trying :) – user565684 Feb 15 at 23:20
• Question has been edited now – user565684 Feb 16 at 23:15

For 1) You can change the code to:

def probably_prime(n, k, a):
"""Return True if n passes k rounds of the Miller-Rabin primality
test (and is probably prime). Return False if n is proved to be
composite.

"""
if n < 2: return False

r, m = 0, n - 1
while m % 2 == 0:
r += 1
m //= 2
for _ in range(k):
#a = randrange(2, n - 1)
x = pow(a, m, n)
if x == 1 or x == n - 1:
continue
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
return True


and then fix some number n, and some k, which is the number of repetitions. Then add at the end:

Result=[]
For a in range(1,n):
Result.append(probably_prime(n,k,a))


Result is a list with all outcomes for $$a$$ in the required range. True means composite (ie witness), and false means liar. To find which $$a$$'s are liars and which are witnesses you can then add:

witnesses=sum(Result)
liars=n-witnesses

• Thank you. I will be working on this in the next few hours. I appreciate your help! – user565684 Feb 15 at 23:31
• Hi Alexandros-I have been trying to get Python code to work, but have not succeeded. When you indicated to fix some number n and k, I am doing as probably_prime(n, k), but you now have probably_prime(n, k, a) so it wants a value for a. How should I be fixing n and k, and isn't a just running through all from 1<a<n-1? Then you have the part to add at the end: Result=[]....I attempted this part after fixing n and k but it did not work, probably because I am not entering the n and k correctly. I you don't mind can you let me know where you think I could be going wrong from what I have indicated. – user565684 Feb 16 at 7:29
• When I said fix n,k i meant to choose a number n which you wish to test whether it is prime or not, and choose a number k as the number of trials. Ie if you wish to test whether n=500 is prime with 5 repetitions, you should add these two lines at the beginning of your code: n=500, k=5. Also, probably_prime(n,k) is not going to work anymore, you need to use probably_prime(n,k,a). – Alexandros Feb 16 at 8:24
• Sorry I got confused as when you indicated to fix n and k you mentioned it after the main part of the code, and that is when I normally enter the definition with the values. So I will now fix at the beginning (by just entering two lines of code n=value,k=value) and add the last two blocks, Result=.. and Witnesses=... at the end of the main part of the code. Thank you again, and I hope you do not mind so many questions :) – user565684 Feb 16 at 15:41
• Alexandros - I still don't have code working. I define n and k at the beginning, then entered the main body of the code, then the last two blocks of code, Lastly, I call it by the definition and it returns "...name a is not defined" Any ideas? I have tried a variety of ways, but just cannot get it working. Maybe it is something simple I am just not aware of. – user565684 Feb 17 at 3:28