## How Fermat Factored Integers

### July 18, 2014

We begin by calculating the table of squares. Instead of printing it as we did on the previous page, we store it in a vector, plus one additional square to act as a sentinal:

`(define limit 10000)`

```(define squares   (let ((squares (make-vector (+ limit 1))))     (do ((i 0 (+ i 1))) ((< limit i) squares)       (vector-set! squares i (* i i)))))```

Then `isqrt`, which calculates the largest number that, when multiplied by itself, doesn’t exceed n, and `in-table`, which returns the square root of n if n is a square or `#f` otherwise, are simple:

```(define (isqrt n)   (do ((i 1 (+ i 1))) ((< n (vector-ref squares i)) (- i 1)))))```

```(define (in-table n)   (let loop ((i 0))     (cond ((< n (vector-ref squares i)) #f)           ((= n (vector-ref squares i)) i)           (else (loop (+ i 1))))))```

We calculate the digital root algorithmically; Fermat would certainly have calculated it visually by casting out nines:

```(define (digital-root n)   (let loop ((n n) (r 0))     (cond ((zero? n) (if (< r 10) r (digital-root r)))           ((< n 10) (loop 0 (+ r n)))           (else (let ((d (modulo n 10)))                   (loop (/ (- n d) 10) (+ r d)))))))```

Here we encode the rules for identifying squares. The `square?` function returns the square root if it exists or `#f` otherwise:

```(define (square? n)   (let* ((ones (modulo n 10))          (tens (modulo (/ (- n ones) 10) 10)))     (and (or (and (= 0 ones) (zero? tens))              (and (= 1 ones) (even? tens))              (and (= 4 ones) (even? tens))              (and (= 5 ones) (= 2 tens))              (and (= 6 ones) (odd? tens))              (and (= 9 ones) (even? tens)))          (member (digital-root n) '(1 4 7 9))          (in-table n))))```

We finish with the function that simulates Fermat’s calculation. First we ensure n fits in our table and remove factors of 2. Then we initialize r and t, iterate until we find a square, and calculate the factors:

```(define (fermat n)   (if (not (< -1 n (* limit limit))) #f     (if (even? n) (list 2 (/ n 2))       (let ((x (isqrt n)))         (if (= (* x x) n) (list x x)           (let loop ((r (- (* x x) n)) (t (+ x x 1)))             (if (not (square? r))                 (loop (+ r t) (+ t 2))                 (let ((x (/ (- t 1) 2)) (y (isqrt r)))                   (list (- x y) (+ x y))))))))))```

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

Pages: 1 2 3