Gaussian Integers, Part 2

November 7, 2014

We assume all the functions of the previous exercise are available. Then the greatest common divisor is computed using this recursive function:

(define (gauss-gcd x y)
(if (gauss-zero? y) x
(let* ((q (gauss-quotient x y))
(r (gauss-remainder x y q)))
(gauss-gcd y r))))

And here are two related functions:

(define (gauss-divides? d n)
(gauss-zero?
(gauss-remainder n d
(gauss-quotient n d))))

(define (gauss-coprime? x y)
(gauss-unit? (gauss-gcd x y)))

Here’s the function that determines if a Gaussian integer is prime. In addition to the three conditions of the description, we ensure the Gaussian integer is not a unit, which is neither prime nor composite:

(define (gauss-prime? x)
(cond ((gauss-unit? x) #f)
((zero? (re x))
(and (prime? (im x))
(= (modulo (im x) 4) 3)))
((zero? (im x))
(and (prime? (re x))
((zero? (im x))
(= (modulo (re x) 4) 3)))
(else (prime? (gauss-norm x)))))

The prime? function may use trial division, or the Miller-Rabin method, or the Baillie-Wagstaff method, or any other to determine if a normal integer is prime. Likewise, when factoring Gaussian integers, any of the normal factorization methods may be applied to the norm. Here’s the function for factoring a Gaussian integer:

(define (gauss-factors x)
(define (find-k p a)
(if (= (expm a (/ (- p 1) 2) p) (- p 1))
(expm a (/ (- p 1) 4) p)
(find-k p (+ a 1))))
(let loop ((g x) (qs (list))
(ps (factors (gauss-norm x))))
(display g) (display qs) (display ps) (newline)
(cond ((null? ps)
(if (gauss-eql? g (gs 1 0)) qs (cons g qs)))
((= (car ps) 2)
(loop (gauss-quotient g (gs 1 1))
(cons (gs 1 1) qs) (cdr ps)))
((= (modulo (car ps) 4) 3)
(let ((q (gs (car ps) 0)))
(loop (gauss-quotient g q)
(cons q qs) (cddr ps))))
(else
(let* ((p (car ps))
(k (find-k p 2))
(q (gauss-gcd (gs p 0) (gs k 1)))
(z (gauss-quotient g q)))
(when (not (gauss-zero? (gauss-remainder g q z)))
(set! q (gauss-conjugate q))
(set! z (gauss-quotient g q)))
(loop z (cons q qs) (cdr ps)))))))

That’s a relatively straight forward implementation of the algorithm described on the previous page. Here’s the example from Jim Lewis:

> (gauss-factors (gauss 361 -1767)
((-7 . -2) (19 . 0) (4 . 1) (2 . 1) (1 . 1))

You can run the program at http://programmingpraxis.codepad.org/ASV1Bh0C, where we used a 2,3,5-wheel for primality checking and integer factorization.