Programming Praxis

To calculate the nth prime we bracket n in loop1, use binary search in loop2, then calculate the prime-counting function and count up or down to the requested n in loop3:

(define (nth-prime n)
  ; (define (approx-pi x) (inexact->exact (round (riemann-r x))))
  (define (approx-pi x) (inexact->exact (round
    (- (logint x) (/ (logint (sqrt x)) 2) (/ (logint (expt x 1/3)) 3)))))
  (if (< n 5) (vector-ref (vector 2 3 5 7) (- n 1))
    (let loop1 ((lo 1) (hi 2))
      (if (< (approx-pi hi) n) (loop1 hi (* hi 2))
        (let loop2 ((lo lo) (hi hi))
          (let* ((mid (quotient (+ lo hi) 2)) (mid-pi (approx-pi mid)))
            (cond ((< n mid-pi) (loop2 lo mid))
                  ((< mid-pi n) (loop2 mid hi))
                  (else (let* ((m (if (prime? mid) mid (prev-prime mid))) (k (prime-pi m)))
                          (cond ((< k n) (let loop3 ((k k) (m m))
                                  (if (= n k) m (loop3 (+ k 1) (next-prime m)))))
                                ((< n k) (let loop3 ((k k) (m m))
                                  (if (= k n) m (loop3 (- k 1) (prev-prime m)))))
                                (else m)))))))))))

The test for n less than 5 is needed because the binary search “steps over” the correct values when n is so small. The call to (prev-prime mid) is necessary because the result of approx-pi is not necessarily prime.

This function makes a very good example of the trade-offs necessary in developing software. The “right” method uses the Riemann function to approximate π, as shown in the commented-out line of code. But that function is called many times in loop1 and loop2, and it’s too slow. A better solution uses the weaker approximation of the first three terms of the Riemann function at the cost of more calls to the next-prime and prev-prime functions but still provides a much faster function.

The millionth prime is 15485863; the billionth prime is 22801763489:

> (nth-prime #e1e6)
15485863
> (nth-prime #e1e9)
22801763489

We used all the machinery of the three previous exercises plus the Baillie-Wagstaff primality checker. You can run the program at http://programmingpraxis.codepad.org/KgVfGvEG.

Pages: 1 2