## Gaussian Integers, Part 1

### November 4, 2014

We represent Gaussian integers as pairs, so *a* + *b i* is represented as `(a . b)`

:

`(define (gauss a b)`

(when (not (integer? a))

(error 'gauss "must be integer"))

(when (not (integer? b))

(error 'gauss "must be integer"))

(gs a b))

The internal function `gs`

does not enforce integer types:

`(define gs cons)`

The pieces of a Gaussian integer can be picked apart with the `re`

and `im`

functions:

`(define (re x) (car x))`

(define (im x) (cdr x))

Scheme offers a complex number type, but doesn’t limit it to integers. Here are conversions to and from Scheme complex numbers:

`(define (gauss-from-complex x)`

(gauss (real-part x) (imag-part x)))

`(define (gauss-to-complex x)`

(make-rectangular (re x) (im x)))

Zero in Gaussian integers is (0 . 0). The role of the unit in regular integers is played by four different values in Gaussian integers, 1, −1, *i* and −*i*. There is no ordering concept in Gaussian integers, but two Gaussian integers can be compared for equality:

`(define (gauss-zero? x)`

(and (zero? (re x)) (zero? (im x))))

`(define (gauss-unit? x)`

(or (and (= (abs (re x)) 1) (zero? (im x)))

(and (zero? (re x)) (= (abs (im x)) 1))))

`(define (gauss-conjugate x)`

(gs (re x) (- (im x))))

`(define (gauss-eql? x y)`

(and (= (re x) (re y))

(= (im x) (im y))))

The *norm* of Gaussian integers corresponds to the absolute value in normal integers:

`(define (gauss-norm x)`

(define (square x) (* x x))

(+ (square (re x)) (square (im x))))

The arithmetic operators follow the descriptions given above:

`(define (gauss-sub . xs)`

(define (sub x y)

(gs (- (re x) (re y)) (- (im x) (im y))))

(cond ((null? xs) (error 'gauss-sub "no operands"))

((null? (cdr xs)) (gauss-negate (car xs)))

(else (let loop ((xs (cdr xs)) (zs (car xs)))

(if (null? xs) zs

(loop (cdr xs) (sub zs (car xs))))))))

`(define (gauss-mul . xs)`

(define (mul x y)

(gs (- (* (re x) (re y))

(* (im x) (im y)))

(+ (* (re x) (im y))

(* (im x) (re y)))))

(let loop ((xs xs) (zs (gs 1 0)))

(if (null? xs) zs

(loop (cdr xs) (mul (car xs) zs)))))

`(define (gauss-quotient num den)`

(let ((n (gauss-norm den))

(r (+ (* (re num) (re den))

(* (im num) (im den))))

(i (- (* (re den) (im num))

(* (re num) (im den)))))

(gs (round (/ r n)) (round (/ i n)))))

`(define (gauss-remainder num den quo)`

(gauss-sub num (gauss-mul den quo)))

You can see these functions in action at http://programmingpraxis.codepad.org/ql2rU4dd.

We’ll discuss Gaussian integers again in the next exercise. In the meantime, you might want to read this description of Gaussian integers by Keith Conrad.

Haskell: http://codepad.org/1VZIecvt

Question1: is the quotient not like this (ax+by) +(ay-bx)i ?

Question2: is it ok to use floor instead of round?

Then in Scala:

class GaussianInteger(val a: Int, val b: Int) {

def addition(that: GaussianInteger): GaussianInteger = new GaussianInteger(a + that.a, b + that.b)

def subtraction(that: GaussianInteger): GaussianInteger = new GaussianInteger(a – that.a, b – that.b)

def crossMultiply(that: GaussianInteger): GaussianInteger = new GaussianInteger(a * that.a – b * that.b, a * that.b + b * that.a)

def quotient(that: GaussianInteger): GaussianInteger = {

val n = that.a * that.a + that.b * that.b

new GaussianInteger((a * that.a – b * that.b) / n, (b * that.a – a * that.b) / n)

}

def remainder(that: GaussianInteger): GaussianInteger = subtraction(quotient(that))

override def toString: String = a + ” + ” + b + “i”

}

The quotient is ((ax+by) + (bx-ay)i) / n; the code was right, I fixed the description. Using floor instead of round doesn’t work, because the norm of the result must be less than half the norm of the divisor, by convention; see the discussion of the Division Theorem in the paper linked from the solution page.

Using Scheme’s complex numbers, I get addition and multiplication for free. I find that Gambit-C does Gaussian rationals exactly, and for it the norm is automatically real. In Guile I had to take the real part explicitly. (I needed the norm to be real in order to test that the remainder is smaller than the divisor.)

(define (zi-conjugate a) (make-rectangular (real-part a) (- (imag-part a))))

(define (zi-norm a) (real-part (* a (zi-conjugate a)))) ;imag-part is zero

(define (zi-quotient a b)

(let ((q (/ (* a (zi-conjugate b)) (zi-norm b))))

(make-rectangular (round (real-part q)) (round (imag-part q)))))

`(define (zi-remainder a b) (- a (* b (zi-quotient a b))))`

Printing the quotient and remainder from the paper:

(write (cons (zi-quotient 27-23i 8+i) (zi-remainder 27-23i 8+i))) (newline)

The Gambit-C interpreter (gsi) prints this, where the absence of decimal points indicates that the components are exact integers (with a zero real component omitted):

(3-3i . -2i)

Guile prints this, where the presence of the decimal points indicates that the computation may have used inexact methods (floating point, I’m sure):

(3.0-3.0i . 0.0-2.0i)

Floating point might fail (rounding, overflowing), but my limited tests, with all components small, were fine even in Guile.