Programming Praxis

We’ve done most of the hard work in previous exercises and the Standard Prelude; see http://programmingpraxis.codepad.org/gjG45ulw for primes, prime?, expm, jacobi, legendre, mod-sqrt, isqrt, divides?, square?, fold-of and list-of. Here is Cornacchia’s algorithm:

(define (cornacchia d p)
  (if (negative? (jacobi (- d) p)) #f
    (let loop ((a p) (b (apply max (mod-sqrt (- d) p))))
      (if (< p (* b b)) (loop b (modulo a b))
        (let* ((n (- p (* b b))) (c (/ n d)))
          (cond ((not (divides? d n)) #f)
                ((not (square? c)) #f)
                (else (list b (isqrt c)))))))))

We can use Cornacchia’s algorithm to solve the two problems given at the beginning of the exercise:

> (cornacchia 4 1733)
(17 19)
> (cornacchia 58 3031201)
(127 228)

To verify Fermat’s Christmas Theorem for the primes less than a thousand, we write a simple list comprehension:

> (list-of (cons p (cornacchia 1 p))
    (p in (primes 1000))
    (= (modulo p 4) 1))
((5 2 1) (13 3 2) (17 4 1) (29 5 2) (37 6 1) (41 5 4)
 (53 7 2) (61 6 5) (73 8 3) (89 8 5) (97 9 4) (101 10 1)
 (109 10 3) (113 8 7) (137 11 4) (149 10 7) (157 11 6)
 (173 13 2) (181 10 9) (193 12 7) (197 14 1) (229 15 2)
 (233 13 8) (241 15 4) (257 16 1) (269 13 10) (277 14 9)
 (281 16 5) (293 17 2) (313 13 12) (317 14 11) (337 16 9)
 (349 18 5) (353 17 8) (373 18 7) (389 17 10) (397 19 6)
 (401 20 1) (409 20 3) (421 15 14) (433 17 12) (449 20 7)
 (457 21 4) (461 19 10) (509 22 5) (521 20 11) (541 21 10)
 (557 19 14) (569 20 13) (577 24 1) (593 23 8) (601 24 5)
 (613 18 17) (617 19 16) (641 25 4) (653 22 13) (661 25 6)
 (673 23 12) (677 26 1) (701 26 5) (709 22 15) (733 27 2)
 (757 26 9) (761 20 19) (769 25 12) (773 22 17) (797 26 11)
 (809 28 5) (821 25 14) (829 27 10) (853 23 18) (857 29 4)
 (877 29 6) (881 25 16) (929 23 20) (937 24 19) (941 29 10)
 (953 28 13) (977 31 4) (997 31 6))

You can run the program at http://programmingpraxis.codepad.org/gjG45ulw.

Pages: 1 2