0

I used some code taken from Rosetta Code. I renamed a few things but I didn't really change anything. 

import random

def is_probable_prime(n, num_trials = 5):
    assert n >= 2
    if n == 2:
        return True
    if n % 2 == 0:
        return False
    s = 0
    d = n-1
    while True:
        quotient, remainder = divmod(d, 2)
        if remainder == 1:
            break
        s += 1
        d = quotient
    assert(2**s * d == n-1)

    def try_composite(a):
        if pow(a, d, n) == 1:
            return False
        for i in range(s):
            if pow(a, 2**i * d, n) == n-1:
                return False
            return True

    for i in range(num_trials):
        a = random.randrange(2, n)
        if try_composite(a):
            return False
    return True

It pretty closely matches some psuedocode. However, when I test the number

123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890123456789012345678901

it returns False. Other (python and java) implementations of Miller-Rabin return True for probable prime. After some testing, try_composite returns True after only 2 rounds! I would really like to know any error, I'm guessing a mis-indent or some feature I don't know about.

qwr
  • 9,525
  • 5
  • 58
  • 102
  • Hi, questions asking us to debug someone's code will generally be closed here. Why don't you try to debug it yourself? Have you ever debugged a program? – sashoalm May 10 '14 at 07:34
  • @sashoalm I have. I see that the first part of the code I believe works as expected - s it 1, d is n/2. When the numbers are passed to try_composite, it seems like neither if statement is triggered, I can't see why because it should work. – qwr May 10 '14 at 07:48
  • I repeat, this site is not intended to provide debugging assistance, so why are you trying to nudge the conversation there? Just debug it, man. Take some time and do your work. See this - https://docs.python.org/2/library/pdb.html - this is the python debugger. – sashoalm May 10 '14 at 07:57

1 Answers1

2

In your try_composite function, the for loop should be for i in range(1,s). Do not test the case where i is zero.

EDIT: Also, you are missing a test in your try_composite function. Here is my version of the pseudocode:

def isPrime(n, k=5):
    def isComposite(s, d):
        x = pow(randrange(2,n-1), d, n)
        if x == 1 or x == n-1: return False
        for r in range(1, s):
            x = pow(x, 2, n)
            if x == 1: return True
            if x == n-1: return False
        return True
    if n < 2: return False
    for p in [2, 3, 5, 7, 11, 13, 17]:
        if n % p == 0: return n == p
    s, d = 0, n-1
    while d % 2 == 0: s, d = s+1, d/2
    for i in range(k):
        if isComposite(s, d): return False
    return True

It's too bad that Python doesn't allow labels on break or continue statements. Here's pseudocode for a much prettier version of the function:

function isPrime(n, k=5)
    if n < 2 then return False
    for p in [2,3,5,7,11,13,17,19,23,29]
        if n % p == 0 then return n == p
    s, d = 0, n-1
    while d % 2 == 0
        s, d = s+1, d/2
    for i from 0 to k
        x = powerMod(randint(2, n-1), d, n)
        if x == 1 or x == n-1 then next i
        for r from 1 to s
            x = (x * x) % n
            if x == 1 then return False
            if x == n-1 then next i
        return False
    return True

Notice the two places where the control flow goes to next i. There is no good way to write that in Python. One choice uses an extra boolean variable that can be set and tested to determine when to bypass the rest of the code. The other choice, which I took above, is to write a local function to perform the task. This "loop-and-a-half" idiom is convenient and useful; it was proposed in PEP 3136 and rejected by Guido.

user448810
  • 17,381
  • 4
  • 34
  • 59