My Weblog

Generating primes in parallel

We can use Miller-Rabin primality testing to test a number is prime or not with certain probability. We can use this method on chunk of data to get the primes in parallel. I wrote this code to see how it goes and it’s quite interesting. My other post about data parallelism and Don Stewart answer regarding different options to improve parallelism. The crux of this code is

```primeRange :: Integer -> Integer -> [ Bool ]
primeRange m n = ( map isPrime [ m .. n ] ) `using` parListChunk 10000 rdeepseq
```

```import Data.Bits
import Control.Parallel.Strategies

powM :: Integer -> Integer -> Integer ->  Integer
powM a d n = powM' a d n where
powM' a d n
| d == 0 = 1
| d == 1 = mod a n
| otherwise = mod q n  where
p = powM'  ( mod ( a^2 ) n ) ( shiftR d 1 ) n
q = if (.&.) d 1 == 1 then mod ( a * p ) n else p

calSd :: Integer -> ( Integer , Integer )
calSd n = ( s , d ) where
s = until ( \x -> testBit ( n - 1) ( fromIntegral x ) )  ( +1 ) 0
d = div ( n - 1 ) (  shiftL 1 ( fromIntegral s )  )

rabinMiller::Integer->Integer->Integer->Integer-> Bool
rabinMiller  n s d a
| n == a = True
| otherwise =  case x == 1 of
True -> True
_ ->   any ( == pred n ) . take ( fromIntegral s )
. iterate (\e -> mod ( e^2 ) n ) \$ x
where
x = powM a d n

isPrime::Integer-> Bool
isPrime n
| n <= 1 = False
| n == 2 = True
| even n = False
| otherwise	= all ( == True ) . map ( rabinMiller n s d ) \$  [ 2 , 3 , 5 , 7 , 11 , 13 , 17 ] where
( s , d ) = calSd n

primeRange :: Integer -> Integer -> [ Bool ]
primeRange m n = ( map isPrime [ m .. n ] ) `using` parListChunk 10000 rdeepseq

main = do
let t = filter ( == True ) . primeRange 2 \$ 10 ^ 6
print . length \$ t

[1 of 1] Compiling Main             ( RabinMiller.hs, RabinMiller.o )
Mukeshs-MacBook-Pro:Haskell mukeshtiwari\$ time ./RabinMiller +RTS -N1
78498

real	0m4.841s
user	0m4.695s
sys	0m0.075s
Mukeshs-MacBook-Pro:Haskell mukeshtiwari\$ time ./RabinMiller +RTS -N2
78498

real	0m3.115s
user	0m5.969s
sys	0m0.215s
Mukeshs-MacBook-Pro:Haskell mukeshtiwari\$ time ./RabinMiller +RTS -N3
78498

real	0m2.687s
user	0m7.556s
sys	0m0.415s
```

I am no expert in parallelism or Haskell so if you have any comment or suggestion then please let me know.

February 5, 2013 -