millerrabin.hs

--------------------------------------------------------------------------------
-- Checks if there is a high probability of the number being prime using the
-- Miller-Rabin primality test
-- modPow must be used with the Integer type so that the multiplication won't
-- overflow because n could be >= sqrt(maxInt)
--------------------------------------------------------------------------------
isPrimeMR :: Int -> Bool
isPrimeMR n
  | n < 2 = False
  | n /= 2 && n `mod` 2 == 0 = False
  | elem n p' = True
  | otherwise = all (test . fromIntegral) $ takeWhile (< n) p'
  where
    p' = drop 1 $ take 8 primes
    n' = fromIntegral n
    n'' = n' - 1
    (u, t) = head $ [(n'' `shiftR` t, t)| t <- [31,30..], n'' `mod` 2 ^ t == 0]
    
    test :: Integer -> Bool
    test a
      = test' 1 x0 (modPow x0 2 n')
      where 
        x0 = modPow a u n'
        test' t' x y
          | t' == t = (y == 1)
          | y == 1 && x /= 1 && x /= n'' = False
          | otherwise = test' (t' + 1) y (modPow y 2 n')