Refactoring a prime number generator generator
Posted: Thu Jan 29, 2009 9:06 pm
I have written a little prime number generator generator (sic), but it is somehow ugly. I think that I need some help refactoring it, or someone pointing out what I could do more concise or efficient.
Description:
The generator returns a function which, whenever called, spits out the next prime number. For example, when it is first called, it returns 2, when it is called again, it returns 3, then 5 etc..
I implemented the sieve of Erastothenes by using a "wheel". This wheel is a list of lists or, to be more precise, it is an alist of numbers associated with their prime factors. Actually, the numbers are associated with a list of the double of each prime factor, since I only look at odd numbers (handling 2 as a special case at the beginning). The wheel always contains all numbers bigger than the current and for them, all prime factors (doubled) smaller than the current.
When the function is called, it first increments current, then checks it against the first number in the wheel. If it is equal, it is not prime, if not, then it is. When it is a prime, I begin to "strike out" all its multiples by putting it into the wheel. When it is not a prime, this part of the wheel is broken up by putting the prime factors on their respective next multiple; then the function recurs, until a prime is found.
Is using an alist sensible here? Should I use a hash table instead?
Description:
The generator returns a function which, whenever called, spits out the next prime number. For example, when it is first called, it returns 2, when it is called again, it returns 3, then 5 etc..
I implemented the sieve of Erastothenes by using a "wheel". This wheel is a list of lists or, to be more precise, it is an alist of numbers associated with their prime factors. Actually, the numbers are associated with a list of the double of each prime factor, since I only look at odd numbers (handling 2 as a special case at the beginning). The wheel always contains all numbers bigger than the current and for them, all prime factors (doubled) smaller than the current.
When the function is called, it first increments current, then checks it against the first number in the wheel. If it is equal, it is not prime, if not, then it is. When it is a prime, I begin to "strike out" all its multiples by putting it into the wheel. When it is not a prime, this part of the wheel is broken up by putting the prime factors on their respective next multiple; then the function recurs, until a prime is found.
Code: Select all
(defun make-prime-generator ()
"returns a function with closure that gives all primes sequentially.
Example: (let ((gimme-prime (make-prime-generator)))
(funcall gimme-prime) ; => 2
(funcall gimme-prime)) ; => 3"
(let ((current 1)
;; wheel is a list of lists. The sublists contain the
;; number and then all prime factors, sorted ascendingly.
(wheel ()))
(labels ((propagate-wheel ()
(do ((wr (cdr wheel))
(wl wheel) ; insert between wl and wr
(p (cdar wheel)))
((null p)
(setf wheel (cdr wheel))
(next-prime))
(let ((strikeout (+ current (car p))))
(cond ((null wr)
(setf (cdr wl) (cons strikeout (car p))
wl (cdr wl)
p (cdr p)))
((< strikeout (caar wr))
(setf (cdr wl) (cons (list strikeout (car p))
wr)
wl (cdr wl)
p (cdr p)))
((= strikeout (caar wr))
(setf (cdar wr) (merge 'list (cdar wr) (list (car p)) #'<))
(setf p (cdr p)))
((> strikeout (caar wr))
(setf wr (cdr wr)
wl (cdr wl)))))))
(put-prime-into-wheel (place number &optional (wl wheel))
(cond ((null wheel) ; outer wheel, when current is 3
(setf wheel (list (list place number))))
((= (caar wl) place) ; inner
(setf (cdar wl) (merge 'list (cdar wl) (list number) #'<)))
((null (cdr wl))
(setf (cdr wl) (list (list place number))))
((< place (caadr wl))
(setf (cdr wl) (cons (list place number) (cdr wl))))
(t
(put-prime-into-wheel place number (cdr wl)))))
(next-prime ()
(case current
(1 (return-from next-prime (incf current)))
(2 (incf current))
(otherwise (incf current 2)))
(if (and (car wheel)
(= (caar wheel) current))
;; not a prime -> rotate wheel
(propagate-wheel) ; will call next-prime again
;; prime -> put it into the wheel, and return it
(progn
(put-prime-into-wheel (* 3 current) (* 2 current))
current))))
#'next-prime)))