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.

Advertisement

Pages: 1 2

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: