Eratosthenes Party Time! (A.k.a. Feedback wanted on this implementation of Eratosthenes Sieve)

So I coudn’t keep my hands off this …

(defn bitset-sieve [n]
  (let [primes (java.util.BitSet. n)
        imax (-> n Math/sqrt int inc)]
    ;; Start by assuming they are all primes
    (.set primes 2 n)
    ;; Now, iterate through the primes.
    (loop [i (.nextSetBit primes 0)]
      ;; Imperative code ...
      (if (or (> i imax)
              (= i -1))
        primes
        (do
          (doseq [j (range (+ i i) n i)]
            ;; Not prime -- multiple of i!
            (.clear primes j))
          (recur (.nextSetBit primes (inc i))))))))

(defn bitset->vec [bs]
  (->> (range (.size bs))
       (filter (fn [x] (.get bs x)))
       (into [])))

With that implementation, I take 20-50 msecs to produce a sieve up to a million. If I mash it into a vector afterwards (bitset->vec), it takes me 1220 msecs, so it seems faster to avoid immutable data structures here. I do start out with all the numbers (range (.size bs)), so that could be sped up, but I think I’d rather avoid creating the vector.

Note that I also use the precomputed prime sequence to jump ahead (recur (.nextSetBit primes (inc i))), not bothering to touch what I know not to be primes.

(defn make-prime-checker [lim]
  (let [sieve (bitset-sieve lim)]
    (fn [n]
      (assert (pos-int? n))
      (assert (< n lim))
      (.get sieve n))))

(time
 (let [checker (make-prime-checker (* 1000 1000))]
   (doseq [i [1001 10321 99097]]
     (prn i 'prime? (checker i))
     )))
;; 1001 prime? false
;; 10321 prime? true
;; 99097 prime? false
;; "Elapsed time: 22.153831 msecs"

Note to self: imperative Clojure code may be fast, and may also hurt the eyes. And I should probably check the result against a list of primes, but now it’s bedtime!

2 Likes