Segmented Sieve Of Eratosthenes
February 5, 2010
We examined the Sieve of Eratosthenes for computing prime numbers less than n in our second exercise nearly a year ago. In today’s exercise, we update that algorithm to compute the prime numbers in a range L to R, where the range doesn’t start at zero and is too large to fit in memory all at once. Specifically, we discuss a function that takes parameters L, for the left end of the range, and R, for the right end of the range, and a parameter B that divides R – L; we also assume that the primes less than the square root of R are known. In practice, L and R will be large, maybe 1016 or 1018, and B will be somewhat smaller, maybe 108 or 1010.
We will explain the algorithm by computing the primes in the range 100 to 200, with B = 10; since we don’t look at even numbers, the algorithm will make five passes, each time sieving twenty numbers.
There are six primes less than the square root of 200: 2, 3, 5, 7, 11, and 13; we will ignore 2, since it is even. Associated with each prime Pk is a number we will call Qk that is the offset within B of the first prime in the current sieving block that is divisible by Pk. For instance, when sieving from 100 to 120, the primes are 3, 5, 7, 11, and 13 and the associated Qs are 2, 2, 2, 10, and 8. Each Qk can be initialized by computing Qk = (-1/2 × (L + 1 + Pk)) mod Pk and re-initialized for the next 2B block by computing Qk = (Qk – B) mod Pk; for instance, the Qs at the second sieving block from 120 to 140 are 1, 2, 6, 0, and 11.
Within each sieving block, we create an array of booleans (in practice, it is common to use bits) that is initialized to true. Then each prime is sieved in the normal way starting at the Qk offset. In our example, prime 3 with associated Q = 2 causes us to reset bits 2, 5, and 8 to false. Likewise, prime 5 causes us to reset bits 2 and 7, prime 7 causes us to reset bits 2 and 9, prime 11 has no odd multiples in the current sieving block, and prime 13 causes us to reset bit 8. After sieving, the bits that are still true are 0, 1, 3, 4, and 6, which correspond to primes 101, 103, 107, 109, and 113.
It is instructive to take five minutes to work out the example sieve by hand.
Your task is to write a function that performs a segmented sieve of Eratosthenes as described above. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
My first try
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 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 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 main = print $ sieve 100 200 10 [3, 5, 7, 11, 13][…] Of Eratosthenes February 5, 2010 Mithrandir Leave a comment Go to comments Today’s Programming Praxis problem describes an algorithm which may be used to find the prime numbers in a large interval, so large […]
Won’t doing “(primes (+ (isqrt r) 1))” sometimes add an extra prime that you won’t need?
For instance, try it with l = 200, and r = 325 (basically, any time r is at least the value of
an even square, but less than an odd square).
I have worked on this problm http://www.spoj.pl/problems/PRIME1/ about this,but this semms too hard for me!!
[…] I would probably choose the first alternative given your problem size. Implementations of both the segmented Sieve of Eratosthenes and Lehmer’s prime-counting function are available at my […]
[…] found a solution here which says the function should be $(-1/2 × (L + 1 + Pk))mod Pk$ , where L:start of range, Pk: […]
[…] can see a simple implementation of a segmented sieve here. Note that a segmented sieve will be very much faster than O’Neill’s priority-queue […]
[…] ver una implementación simple de un tamiz segmentado aquí. Tenga en cuenta que un tamiz segmentado será mucho más rápido que el tamiz de cola de prioridad […]
[…] can see a simple implementation of a segmented sieve here. Note that a segmented sieve will be very much faster than O’Neill’s priority-queue […]