MillerRabin.hs: Probabilistic Primality Tests in Haskell

I've been getting re-acquainted with Haskell, so I decided to try my hand at implementing Miller-Rabin primality test¹. The functions millerRabin and millerRabinIO are included in the module MillerRabin, which can be downloaded here. A fast modular exponentiation function powMod is also included.

The Miller-Rabin primality test is a way of quickly determining whether a given number $n$ is probably prime or composite. Each prime is checked against $k$ randomly chosen integers $b_1 … b_k$, each of which can be used to show that $n$ is composite, but not that it is prime. If $b_i$ is used to show that $n$ is composite, we call it a witness to the compositeness of $n$. If If $n$ is composite, then the probability that a randomly chosen number $b_i$ is not a witness has an upper bound of $\frac{1}{4}$, so if we check $b_1$ through $b_k$, the probability that we will find no witnesses is less than $\frac{1}{4^k}$.

I have written the module to be consistent with the style of Haskell's System.Random. As with most functions in System.Random, there are two versions of the Miller-Rabin test:

millerRabin :: (Integral a, Random a, RandomGen g) 
            => a -> Int -> g -> (Bool, g)

millerRabinIO :: (Integral a, Random a) => a -> Int -> IO Bool

millerRabin enables the programmer to provide their own random number generator, while millerRabinIO uses the global random number generator.

Example usage:

λ> :load MillerRabin
[1 of 1] Compiling MillerRabin      ( MillerRabin.hs, interpreted )
Ok, one module loaded.
λ> millerRabinIO (2^127 - 1) 4
True
λ> millerRabinIO (2^257 - 1) 4
False
λ> import System.Random
λ> millerRabin 188748146801 4 (mkStdGen 1729)
(True,1689113473 525453832)
λ> millerRabin 443372888629441 4 (mkStdGen 1729)
(False,1689113473 525453832)

These functions are intended to be used as strategies in a larger isPrime function. They do not save time by checking if $n$ is even or is a small prime.

The Miller-Rabin primality test can be modified to make it deterministic for numbers smaller than some threshold; instead of choosing $b_i$ randomly, $b_1 … b_k$ is a pre-determined set known to catch all composite numbers less than some bound. OEIS A014233 is the smallest odd number for which the test fails to reveal compositeness when checked against the first $k$ primes; this sequence shows that it is possible to check compositeness for all numbers up to 2⁶⁴ by checking against the first 13 primes. As a further improvement on that, the website miller-rabin.appspot.com shows the smallest sets of bases that have been found so far. I haven't yet included a deterministic Miller-Rabin function in MillerRabin.hs because I haven't figured out to make it polymorphic over both Ints and Integers; supposedly the AKS primality test is more commonly used when a deterministic primality test is needed, so I may implement that instead.