Today’s Programming Praxis problem describes an algorithm which may be used to find the prime numbers in a large interval, so large that it cannot be kept in memory at once. See their description there, I’ll only give the code in Haskell.

First, some utility functions:

`firstQ :: Int -> Int -> Int`

firstQ l pk = (-(1+l+pk) `div` 2) `mod` pk

nextQ :: Int -> Int -> Int -> Int

nextQ b pk qk = (qk - b) `mod` pk

Then we sieve one block and convert the values to the primes:

`sieveBlock :: [Int] -> [Int] -> [Int] -> [Int]`

sieveBlock [] [] blk = blk

sieveBlock (p:ps) (o:ofs) blk = filter f (sieveBlock ps ofs blk)

where

f x = (x < o) || ((x - o) `mod` p /= 0)

getValsfromBlock :: Int -> [Int] -> [Int]

getValsfromBlock l bl = map (\x->l+1+2*x) bl

Now, it’s time to obtain all the numbers. This post will make use of the until function which is very similar to what some will use for the `for` construct in C.

`sieve :: Int -> Int -> Int -> [Int] -> [Int]`

sieve l r b primes = middle $ until stop incr start

where

middle (_, v, _) = v

stop (sb, _, _) = sb == r

start = (l, [], map (firstQ l) primes)

incr (startblock, primessofar, q) = (nextblock, nextprimes, nextq)

where

array = [0 .. b-1]

nextblock = startblock + 2 * b

nextprimes = primessofar ++ (getValsfromBlock startblock bl)

nextq = zipWith (nextQ b) primes q

bl = sieveBlock primes q array

The main is very simple after we know about the print function and that it is a monadic one

`main = print $ sieve 100 20000 10 [3, 5, 7, 11, 13]`

That’s all for now.