### Everybody loves the Sieve of Eratosthenes

unsorted — cgrand, 30 July 2009 @ 11 h 27 min

If I judge by traffic logs for this blog, many newcomers want to compute prime numbers in Clojure.

A recent thread on the mailing list prompted me to test an idea I had for a while: to use a map to implement the Sieve of Eratosthenes.

```(defn primes [max]
(let [enqueue (fn [sieve n factor]
(let [m (+ n factor)]
(assoc sieve m
(conj (sieve m) factor))))
next-sieve (fn [sieve candidate]
(if-let [factors (sieve candidate)]
(reduce #(enqueue %1 candidate %2)
(dissoc sieve candidate)
factors)
(enqueue sieve candidate candidate)))]
(apply concat (vals (reduce next-sieve {} (range 2 max))))))
```

where the sieve is a map from the next non-prime numbers to their factors. It’s so naive that even numbers are tested for primality but it doesn’t perform that bad: on my box, it takes 3s to compute all primes below 1,000,000 (it’s roughly as fast as clojure.contrib.lazy-seqs/primes when the seq isn’t yet memoized).

I wasn’t happy with the way I handled the case where a non-prime was already checked off (ie was a key of the map): I was conjing onto the list of prime factors for this number. If instead I tried to check off n+p (where n is the already known non-prime and p the current prime) or n+2p or n+3p… until I found a yet unknown non-prime I wouldn’t need to maintain a list of factors, nor to conj or reduce over it. And I was hoping that less allocations would yield better perfs.

Here is the second iteration:

```(defn primes2 [max]
(let [enqueue (fn [sieve n factor]
(let [m (+ n factor)]
(if (sieve m)
(recur sieve m factor)
(assoc sieve m factor))))
next-sieve (fn [sieve candidate]
(if-let [factor (sieve candidate)]
(-> sieve
(dissoc candidate)
(enqueue candidate factor))
(enqueue sieve candidate candidate)))]
(vals (reduce next-sieve {} (range 2 max)))))
```

and it computes all the primes below 1,000,000 in 1.8s instead of 3s but it still tests even numbers for primality.

primes3 is primes2 modified to only test odd numbers:

```(defn primes3 [max]
(let [enqueue (fn [sieve n factor]
(let [m (+ n (+ factor factor))]
(if (sieve m)
(recur sieve m factor)
(assoc sieve m factor))))
next-sieve (fn [sieve candidate]
(if-let [factor (sieve candidate)]
(-> sieve
(dissoc candidate)
(enqueue candidate factor))
(enqueue sieve candidate candidate)))]
(cons 2 (vals (reduce next-sieve {} (range 3 max 2))))))
```

and it computes the same list of primes in 1.5s.

Out of curiosity, I wrote a lazy version of primes3:

```(defn lazy-primes3 []
(letfn [(enqueue [sieve n step]
(let [m (+ n step)]
(if (sieve m)
(recur sieve m step)
(assoc sieve m step))))
(next-sieve [sieve candidate]
(if-let [step (sieve candidate)]
(-> sieve
(dissoc candidate)
(enqueue candidate step))
(enqueue sieve candidate (+ candidate candidate))))
(next-primes [sieve candidate]
(if (sieve candidate)
(recur (next-sieve sieve candidate) (+ candidate 2))
(cons candidate
(lazy-seq (next-primes (next-sieve sieve candidate)
(+ candidate 2))))))]
(cons 2 (lazy-seq (next-primes {} 3)))))
```

and, surprisingly (better locality?), it computes the primes below 1,000,000 in 1s.

1. Here's my naive, ugly and untested take on the problem:
(defn not-divisible-by?[num denum]
(not (= (mod num denum) 0)))

(defn primes
([pos] (primes (vec (range 2 pos)) 0))
([clist pos]
(cond
(>= pos (count clist)) clist
:else
(let [val (nth clist pos)
starting (inc pos)]
(recur (into (subvec clist 0 starting) (vec (filter #(not-divisible-by? % val) (subvec clist starting (count clist))))) starting)))))

Comment by zemariamm — 15 September 2009 @ 4 h 36 min
2. [...] reading this post [...]

3. [...] favorite solution was by Cristophe Grande, not just because he followed a similar approach to me, but also because he was very clever in his [...]

4. [...] of SoE. I would recommend checking out Christophe Grand’s treatise on the subject titled Everybody loves the Sieve of Eratosthenes for a great discussion on writing real world prime sieves in [...]

5. [...] Clojure is very big on laziness… and I’m feeling to lazy to re-implement the Sieve in a one-liner. I’m sure it could be done! In the meantime, I’ll point you to cgrand’s great article about the Sieve of Eratosthenes in clojure. [...]

6. thanks for this good article..

Comment by heiji — 27 August 2011 @ 16 h 15 min
7. [...] I wasn’t happy with this performance, especially after seeing Christophe Grand’s solution and timing on his blog post Everybody loves the Sieve of Eratosthenes. [...]

8. アバクロ 通販

Comment by アバクロ パンツ — 11 September 2013 @ 12 h 06 min
9. […] You should follow the link. There are many examples and answers there. Also clj-me.cgrand.net/2009/07/30/…. […]

10. My solution cheats a little in that it doesn’t compute the next prime number that will be used to make multiples. Instead it uses all odd integers between 2 and the limit.

(defn primes [below]
(remove (set (mapcat #(range (* % %) below %)
(range 3 (Math/sqrt below) 2)))
(cons 2 (range 3 below 2))))

Comment by Thomas Varney — 2 December 2015 @ 5 h 01 min