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.