I just got into haskell.
I'm trying to write a prime number generator that would return me a prime number depending on the number of bits I give it. Do I use some sort of Sieve of Eratosthenes? Would this be the quickest way? Currently, I already have a prime number checker for Miller-Rabin. Is there a correct way and a wrong way to do this? Also, I want to be able to generate large numbers very quickly.
Ex. generate a prime number that is of 32 bits
genp 32
Code so far.
import System.IO
import System.Random
import Data.List
import Control.Monad
import Control.Arrow
primesTo100 = [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97]
powerMod :: (Integral a, Integral b) => a -> a -> b -> a
powerMod m _ 0 = 1
powerMod m x n | n > 0 = join (flip f (n - 1)) x `rem` m where
f _ 0 y = y
f a d y = g a d where
g b i | even i = g (b*b `rem` m) (i `quot` 2)
| otherwise = f b (i-1) (b*y `rem` m)
witns :: (Num a, Ord a, Random a) => Int -> a -> IO [a]
witns x y = do
g <- newStdGen
let r = [9080191, 4759123141, 2152302898747, 3474749600383, 341550071728321]
fs = [[31,73],[2,7,61],[2,3,5,7,11],[2,3,5,7,11,13],[2,3,5,7,11,13,17]]
if y >= 341550071728321
then return $ take x $ randomRs (2,y-1) g
else return $ snd.head.dropWhile ((<= y).fst) $ zip r fs
isMillerRabinPrime :: Integer -> IO Bool
isMillerRabinPrime n | n `elem` primesTo100 = return True
| otherwise = do
let pn = pred n
e = uncurry (++) . second(take 1) . span even . iterate (`div` 2) $ pn
try = return . all (\a -> let c = map (powerMod n a) e in
pn `elem` c || last c == 1)
witns 100 n >>= try
type Prime = Integer
isProbablyPrime :: Prime -> Bool
isProbablyPrime n = isMillerRabinPrime n
pickFirstFrom :: Integer -> Prime
pickFirstFrom n = head $ filter isProbablyPrime [n..]
numBits = 1024
constantStdGen = mkStdGen 1234567 -- a random number
randomByBits n = fst $ randomR (2^(n-1), 2^n-1) constantStdGen
answer = pickFirstFrom (randomByBits numBits)
Using Library functions with RSA:
import Control.Monad.Fix
import Math.NumberTheory.Primes
rndPrime :: Int -> IO Integer
rndPrime bits =
fix $ \again -> do
x <- fmap (.|. 1) $ randomRIO (2^(bits - 1), 2^bits - 1)
if isPrime x then return x else again
rndPrimes :: Int -> IO (Integer, Integer)
rndPrimes bits = do
p <- rndPrime bits
fix $ \again -> do
q <- rndPrime bits
if p /= q then return (p, q) else again
Thanks guys, I really appreciate your help.