Miller Rabin Primality test

Question

I have wrote a Miller Rabin primality test in C sharp, but it return false on every input.

Here's the code:

    static Boolean MilRab(UInt64 n)
    {
        UInt64[] ar = new UInt64[] { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29 };
        for (int i = 0; i < 10; i++)
        {
            if (Tanu(ar[i], n) == true) return false;
        }
        return true;
    }//MilRab

    static Boolean Tanu(UInt64 a, UInt64 n)
    {
        UInt64 b = n - 1;
        UInt64 d = 1;
        UInt64 x;
        Int16 i;
        for (i = 63; i >= 0; i--) if (((b >> i) & 1) == 1) break;

        for (;i>=0;i--)
        {
            x = d;
            d = ((d * d) % n);
            if (d == 1 && x != 1 && x != n - 1) return true;
            if (b>>i == 1) d = (d * a) % n;
        }
        if (d != 1) return true;
      return false;
    }//Tanu

What do you think the problem can be? I spent a hole day debugging and its driving me crazy. Thank you.

Answer 1

As Amir said,

http://rosettacode.org/wiki/Miller-Rabin_primality_test#C.23

has the solution. But.. I tried it and even for very small input values like 79, it gives wrong results (79 is prime but comes out repeatedly as not prime). The reason is an overflow in the power (e. g. 11^78 is something like 1E+81, which is not representable by a 32 or 64bit integer).

So here

public static class RabinMiller
{
    public static bool IsPrime(int n, int k)
    {
        if(n < 2)
        {
            return false;
        }
        if(n != 2 && n % 2 == 0)
        {
            return false;
        }
        int s = n - 1;
        while(s % 2 == 0)
        {
            s >>= 1;
        }
        Random r = new Random();
        for (int i = 0; i < k; i++)
        {
            double a = r.Next((int)n - 1) + 1;
            int temp = s;
            int mod = (int)Math.Pow(a, (double)temp) % n;
            while(temp != n - 1 && mod != 1 && mod != n - 1)
            {
                mod = (mod * mod) % n;
                temp = temp * 2;
            }
            if(mod != n - 1 && temp % 2 == 0)
            {
                return false;
            }
        }
        return true;
    }
}

you should add a function

    static int ModuloPower(int a, int b, int n)
    {
        // return (a^b)%n
        int res = 1;
        for (int i = 0; i < b; ++i)
            res = (res * a) % n;
        return res;
    }

and change the power/modulo line to

    int mod = ModuloPower(a, temp, n); // (int)Math.Pow(a, (double)temp) % n;

This can of course be optimized again.. see code examples for other languages on the rosettacode site for some more clever ways to calculate the power/modulo faster.

(I tried editing code there but would have to sign up, so I post here.)

Answer 2

Check this implementation for int and BigInteger

http://rosettacode.org/wiki/Miller-Rabin_primality_test#C.23