Miller-Rabin Recursive Primality Test

Question

I'm working on a primality test and have written a recursive function that returns the value of the function

\$b^{q-1} \bmod q\$

where \$3<= q <= 32000\$

Is there any way to speed up my function? It works, but takes a while to return the answer as \$q\$ approaches 32000.

Variables:

pow = \$q-1\$

mod = \$q\$

b is a variable ranging from \$1 < b < q \$

If q is prime, then b will be = to q if not, b will be a "strong" feature of non-primality. See Miller–Rabin primality test.

public int fF(int q)
    {
        int b = 2, v = 0;
        while(b < q)
        {
            v = operate(b, q-1, q);
            if (v != 1)
                break;
            b++;
        }
        return b;
    }

 int operate(int b, int pow, int mod)
    {
        if (pow == 2)
           return (b * b) % mod;
        return (pow % 2 != 0) ? (b * operate(b, pow - 1, mod)) % mod : (operate(b, pow / 2, mod) * operate(b, pow / 2, mod)) % mod;
    }

Answer 1

Naming

Oh my, what would Mr.Maintainer think if he would inherit the code... single letter variable names, methodnames like fF. He would have a hard time to figure out what is happening.

So let us clean this a little bit

public int fF(int possiblePrime)
{
    int baseNumber = 2, v = 0;

    while (baseNumber < possiblePrime)
    {
        int exponent = possiblePrime - 1;
        v = operate(baseNumber, exponent, possiblePrime);
        if (v != 1)
            break;
        baseNumber++;
    }
    return baseNumber;
}

int operate(int baseNumber, int exponent, int divisor)
{
    if (exponent == 2)
        return (baseNumber * baseNumber) % divisor;
    return (exponent % 2 != 0) ? (baseNumber * operate(baseNumber, exponent - 1, divisor)) % divisor : (operate(baseNumber, exponent / 2, divisor) * operate(baseNumber, exponent / 2, divisor)) % divisor;
}

Style

As many will agree, using braces {}, for single if statements also, is a have to.So let us use them and also let us remove the tenary expression and add an int result which we will return

public int fF(int possiblePrime)
{
    int baseNumber = 2, v = 0;

    while (baseNumber < possiblePrime)
    {
        int exponent = possiblePrime - 1;
        v = operate(baseNumber, exponent, possiblePrime);
        if (v != 1)
        {
            break;
        }
        baseNumber++;
    }
    return baseNumber;
}

int operate(int baseNumber, int exponent, int divisor)
{
    int result = 0;
    if (exponent == 2)
    {
        result = (baseNumber * baseNumber) % divisor;
    }
    else if (exponent % 2 != 0)
    {
        result = (baseNumber * operate(baseNumber, exponent - 1, divisor)) % divisor;
    }
    else
    {
        result = (operate(baseNumber, exponent / 2, divisor) * operate(baseNumber, exponent / 2, divisor)) % divisor;
    }
    return result;
}

Refactoring

Now let us focus on operate()

What you are doing, is always calling number * number % divisor so let us extract this to a method

private int calculateProductModulo(int firstValue, int secondValue, int moduloNumber)
{
    return (firstValue * secondValue) % moduloNumber;
}  

The operate() method now looks

int operate(int baseNumber, int exponent, int divisor)
{
    int result = 0;
    if (exponent == 2)
    {
        result = calculateProductModulo(baseNumber, baseNumber, divisor);
    }
    else if (exponent % 2 != 0)
    {
        result = calculateProductModulo(baseNumber, operate(baseNumber, exponent - 1, divisor), divisor);
    }
    else
    {
        result = calculateProductModulo(operate(baseNumber, exponent / 2, divisor), operate(baseNumber, exponent / 2, divisor), divisor);
    }
    return result;
}   

If we now extract the recursive calls out of the call to calculateProductModulo() we will see clearly what you have stated in your answer

int operate(int baseNumber, int exponent, int divisor)
{
    int result = 0;
    if (exponent == 2)
    {
        result = calculateProductModulo(baseNumber, baseNumber, divisor);
    }
    else if (exponent % 2 != 0)
    {
        int recursiveResult = operate(baseNumber, exponent - 1, divisor);
        result = calculateProductModulo(baseNumber, recursiveResult, divisor);
    }
    else
    {
        int recursiveResult1 = operate(baseNumber, exponent / 2, divisor);
        int recursiveResult2 = operate(baseNumber, exponent / 2, divisor);
        result = calculateProductModulo(recursiveResult1, recursiveResult2, divisor);
    }
    return result;
}   

The code is calling 2 times the same method with the same arguements.

Let us eleminate the double calling

int operate(int baseNumber, int exponent, int divisor)
{
    int result = 0;
    int recursiveResult = 0
    if (exponent == 2)
    {
        result = calculateProductModulo(baseNumber, baseNumber, divisor);
    }
    else if (exponent % 2 != 0)
    {
        recursiveResult = operate(baseNumber, exponent - 1, divisor);
        result = calculateProductModulo(baseNumber, recursiveResult, divisor);
    }
    else
    {
        recursiveResult = operate(baseNumber, exponent / 2, divisor);
        result = calculateProductModulo(recursiveResult , recursiveResult , divisor);
    }
    return result;
}   

Answer 2

I've realized my problem, so I've changed

return (pow % 2 != 0) ? (b * operate(b, pow - 1, mod)) % mod : (operate(b, pow / 2, mod) * operate(b, pow / 2, mod)) % mod; 

to

return (pow % 2 != 0) ? (b * operate(b, pow - 1, mod)) % mod : (int)Math.Pow(operate(b, pow / 2, mod)),2) % mod;

In the first return I was calling the recursive function twice and then evaluating the square. Instead, in the second, I call it once and then evaluate the square. It runs much faster now.

Answer 3

The idea of fast exponentiation is corrupted by the following statement

result=operate(b, pow / 2, mod) * operate(b, pow / 2, mod)       #  (1)

instead of

result=operate(b, pow / 2, mod)**2

or

aux=operate(b, pow / 2, mod)
result=(aux*aux)%mod

It actually slows down the performance from \$O(log(\text{pow})\$ multiplication to \$pow-1\$ multiplications. This is the performance of the dumb exponentiation algorithm (multiplying \$b\$ \$e-1$ times by itself. The performance gain comes from avoiding this second evaluation in (1).