2012-03-22

SICP Exercise 2.86: Complex Complex Numbers

Suppose we want to handle complex numbers whose real parts, imaginary parts, magnitudes, and angles can be either ordinary numbers, rational numbers, or other numbers we might wish to add to the system. Describe and implement the changes to the system needed to accommodate this. You will have to define operations such as sine and cosine that are generic over ordinary numbers and rational numbers.

Type Validation

If you recall, back in exercise 2.83 we augmented all of the procedures installed under the 'make key to validate the (Scheme built-in) types of the numeric operands passed to them to ensure that we could only construct valid representations of numbers. In the case of the complex package we restricted the real and imaginary parts, magnitudes and angles to Scheme real numbers by using real? to test each operand passed to the make-from-real-imag and make-from-mag-ang procedures defined in the rectangular and polar packages.

I'm assuming that, as we're going to allow the use of the integer, rational and real number representations provided by the integer, rational and real packages as the numeric operands for constructing complex representations, we're now going to disallow the use of the built-in Scheme real numbers. This means we need to replace the existing tests with ones that test that the operands are using numeric representations from within our tower-of-types. To ensure that we're not trying to build complex numbers with components that are themselves complex numbers we'll go one step further and ensure that the type we're using is lower in the tower than 'complex.

In the preceding exercise we devised a simple test for the first part of this restriction, which we used in apply-generic to determine whether or not we should drop the result down the tower-of-types: we check that the value is a pair and, if it is, get its type-tag and check that this is in the tower-of-types using memq. We can extract this to its own procedure:
(define (in-tower? value)
  (and (pair? value) (memq (type-tag value) tower-of-types)))
To perform the second part of the test, i.e. to check that the value's type is lower than 'complex, we can use the fact that the tower-of-types is a list ordered from lowest to highest in the tower. Recall that if the procedure memq finds the requested element in the provided list then it returns the tail of that list from that element onwards. If we use memq with the type we want to ensure is higher in the tower ('complex in this case, but let's call it T) then, provided T is present in the tower-of-types, memq will return a list consisting of T and all higher types. We can then check whether or not the type of the value we want to check is lower than T by testing to see whether it's present in the list of types that are T or higher. Iff it's not present in that list then it's a lower type.

Here's our procedure:
(define (is-lower? value type)
  (let ((type-and-higher (memq type tower-of-types)))
    (if (and type-and-higher
             (in-tower? value))
        (not (memq (type-tag value) type-and-higher))
        (error "Either value's type or type is not in tower-of-types"
               (list value type)))))
Let's see it in action:
> (is-lower? (make-integer 4) 'complex)
#t
> (is-lower? (make-rational 3 4) 'real)
#t
> (is-lower? (make-rational 3 4) 'rational)
#f
> (is-lower? (make-rational 3 4) 'integer)
#f
We can then update the rectangular and polar packages to use this test instead of real?:
(define (install-rectangular-package)
  ;; internal procedures
  …
  (define (make-from-real-imag x y)
    (if (and (is-lower? x 'complex) (is-lower? y 'complex))
        (cons x y)
        (error "non-real real or imaginary value" (list x y))))
  …
  (define (make-from-mag-ang r a) 
    (if (and (is-lower? r 'complex) (is-lower? a 'complex))
        (cons (* r (cos a)) (* r (sin a)))
        (error "non-real magnitude or angle" (list r a))))
  …
  ;; interface to the rest of the system
  …
  'done)

(define (install-polar-package)
  ;; internal procedures
  …
  (define (make-from-mag-ang r a)
    (if (and (is-lower? r 'complex) (is-lower? a 'complex))
        (cons r a)
        (error "non-real magnitude or angle" (list r a))))
  …
  (define (make-from-real-imag x y) 
    (if (and (is-lower? x 'complex) (is-lower? y 'complex))
        (cons (sqrt (add (square x) (square y)))
              (atan y x))
        (error "non-real real or imaginary value" (list x y))))
  …
  ;; interface to the rest of the system
  …
  'done)

Generic Arithmetic Operations

Next let's think about the further changes we need to make to the rectangular and polar packages, as well as the changes we'll need to make to the complex package, in order that the operations defined within these packages can work with our system's numeric representations.

The exercise says "you will have to define operations such as sine and cosine that are generic over ordinary numbers and rational numbers." We'll need to go further than just defining generic equivalents for the trigonometric functions though. We'll need to define generic equivalents for all mathematical functions used within the complex packages. Once we've done that we can then replace all of the built-in mathematical functions within the complex packages with their generic equivalents.

Looking at the complex package and the two sub-packages, rectangular and polar we can see that the built-in mathematical functions that we'll need generic equivalents for are: +, -, *, /, =, sqrt, square, atan, cos and sin. We've already got generic equivalents for the first five: add, sub, mul and divide respectively for the first four, and equ?/=zero? for = depending upon the value we're comparing.

Let's move on to the remaining functions.

We can produce the generic operations in the usual manner:
(define (square-root x) (apply-generic 'sqrt x))
(define (square x) (apply-generic 'square x))
(define (arctan x y) (apply-generic 'arctan x y))
(define (cosine x) (apply-generic 'cosine x))
(define (sine x) (apply-generic 'sine x))
All that remains is to install appropriate implementations of these in each of the packages lower in the tower-of-types than 'complex.

The procedures we need to install for 'square are straightforward. For 'integer and 'real values we can just use the built-in * operation and multiply the underlying numeric value by itself. For 'rational numbers we can simply use mul-rat to multiply the rational number by itself. The result can then be tagged as normal.

As for the square root and trigonometric functions, while we could calculate these using the techniques we learned back in section 1.3.3 and its associated exercises (e.g. such as exercise 1.39), we'll go easy on ourselves and just delegate to the corresponding built-in operations appropriately for each type. However, we do need to treat the results of these operations carefully:
  • Given Scheme integer, rational or real values, the results of the built-in trigonometric operations may themselves be Scheme real numbers so, apart from in the real package, we can't just tag the results of the delegated operation normally. Instead we'll have to create real representations from them.
  • Worse still, the result of the built-in sqrt operation may be a Scheme complex number. In this case we'll need to extract out the real and imaginary parts of the result and use these to create a complex representation... And as we're updating the system so that complex numbers can only be constructed using our system's own number representations, we'll have to turn the real and imaginary parts into real representations in order for us to construct the complex representation.
Once we've constructed the appropriate representation for the results we can then rely upon the work we did in exercise 2.85 to drop the number's type to the simplest representation available.

It's possibly worth noting that the complexity of dealing with the results of sqrt points to another piece of work we need to do. As complex representations can only be created using representations from within the tower-of-types we'll need to update the coercion procedure real->complex to that it too turns the components into real representations.

Given all that, here are the changes required to the integer, rational and real packages in order to install the procedures and support the complex package restrictions:
(define (install-integer-package)
  …
  (put 'sqrt '(integer)
       (lambda (x)
         (let ((root (sqrt x)))
           (make-complex-from-real-imag (make-real (real-part root))
                                        (make-real (imag-part root))))))
  (put 'square '(integer)
       (lambda (x) (tag (* x x))))
  (put 'arctan '(integer integer)
       (lambda (x y) (make-real (atan x y))))
  (put 'cosine '(integer)
       (lambda (x) (make-real (cos x))))
  (put 'sine '(integer)
       (lambda (x) (make-real (sin x))))
  …
  'done)

(define (install-rational-package)
  ;; internal procedures
  …
  (define (sqrt-rat x)
    (let ((root (sqrt (/ (numer x) (denom x)))))
      (make-complex-from-real-imag (make-real (real-part root))
                                   (make-real (imag-part root)))))
  (define (square-rat x)
    (mul-rat x x))
  (define (arctan-rat x y)
    (atan (/ (numer x) (denom x))
          (/ (numer y) (denom y))))
  (define (cosine-rat x)
    (cos (/ (numer x) (denom x))))
  (define (sine-rat x)
    (sin (/ (numer x) (denom x))))
  …
  ;; interface to rest of the system
  …
  (put 'sqrt '(rational)
       (lambda (x) (make-real (sqrt-rat x))))
  (put 'square '(rational)
       (lambda (x) (tag (square-rat x))))
  (put 'arctan '(rational rational)
       (lambda (x y) (make-real (arctan-rat x y))))
  (put 'cosine '(rational)
       (lambda (x) (make-real (cosine-rat x))))
  (put 'sine '(rational)
       (lambda (x) (make-real (sine-rat x))))
  …
  'done)

(define (install-real-package)
  …
  (define (real->complex r) (make-complex-from-real-imag (tag r) (tag 0)))
  …
  (put 'sqrt '(real)
       (lambda (x)
         (let ((root (sqrt x)))
           (make-complex-from-real-imag (tag (real-part root))
                                        (tag (imag-part root))))))
  (put 'square '(real)
       (lambda (x) (tag (* x x))))
  (put 'arctan '(real real)
       (lambda (x y) (tag (atan x y))))
  (put 'cosine '(real)
       (lambda (x) (tag (cos x))))
  (put 'sine '(real)
       (lambda (x) (tag (sin x))))
  …
  'done)
Before we go any further, let's give this a spin:
> (square-root (make-real 16.0))
(integer . 4)
> (square (make-rational 3 5))
(rational 9 . 25)
> (cosine (make-integer 0))
(integer . 1)
> (cos 0)
1.0
> (sine (make-integer 3))
(rational 5084384125703515 . 36028797018963968)
> (/ 5084384125703515.0 36028797018963968.0)
0.141120008059867
> (sin 3)
0.141120008059867
> (arctan (make-integer 3) (make-integer 4))
(rational 5796142707547873 . 9007199254740992)
> (/ 5796142707547873.0 9007199254740992.0)
0.643501108793284
> (atan 3 4)
0.643501108793284

Updating the complex Packages

All that remains is to go through the rectangular, polar and complex packages and replace all the built-in operations with our generic equivalents.

Note that the changes we'll make to make-from-mag-ang in the rectangular package and to make-from-real-imag in the polar package will have the effect of dropping the calculated components of the created complex representation. To keep everything consistent, let's explicitly drop the components passed to make-from-real-imag in the rectangular package and to make-from-mag-ang in the polar package.

Here are the changes to the packages:
(define (install-rectangular-package)
  ;; internal procedures
  …
  (define (make-from-real-imag x y)
    (if (and (is-lower? x 'complex) (is-lower? y 'complex))
        (cons (drop x) (drop y))
        (error "non-real real or imaginary value" (list x y))))
  (define (magnitude z)
    (square-root (add (square (real-part z))
                 (square (imag-part z)))))
  (define (angle z)
    (arctan (imag-part z) (real-part z)))
  (define (make-from-mag-ang r a) 
    (if (and (is-lower? r 'complex) (is-lower? a 'complex))
        (cons (mul r (cosine a)) (mul r (sine a)))
        (error "non-real magnitude or angle" (list r a))))
  …
  ;; interface to the rest of the system
  …
  'done)

(define (install-polar-package)
  ;; internal procedures
  …
  (define (make-from-mag-ang r a)
    (if (and (is-lower? r 'complex) (is-lower? a 'complex))
        (cons (drop r) (drop a))
        (error "non-real magnitude or angle" (list r a))))
  (define (real-part z)
    (mul (magnitude z) (cosine (angle z))))
  (define (imag-part z)
    (mul (magnitude z) (sine (angle z))))
  (define (make-from-real-imag x y) 
    (if (and (is-lower? x 'complex) (is-lower? y 'complex))
        (cons (square-root (add (square x) (square y)))
              (arctan y x))
        (error "non-real real or imaginary value" (list x y))))
  …
  ;; interface to the rest of the system
  …
  'done)

(define (install-complex-package)
  ;; imported procedures from rectangular and polar packages
  …
  ;; internal procedures
  …
  (define (add-complex z1 z2)
    (make-from-real-imag (add (complex-real-part z1) (complex-real-part z2))
                         (add (complex-imag-part z1) (complex-imag-part z2))))
  (define (sub-complex z1 z2)
    (make-from-real-imag (sub (complex-real-part z1) (complex-real-part z2))
                         (sub (complex-imag-part z1) (complex-imag-part z2))))
  (define (mul-complex z1 z2)
    (make-from-mag-ang (mul (complex-magnitude z1) (complex-magnitude z2))
                       (add (complex-angle z1) (complex-angle z2))))
  (define (div-complex z1 z2)
    (make-from-mag-ang (div (complex-magnitude z1) (complex-magnitude z2))
                       (sub (complex-angle z1) (complex-angle z2))))
  (define (equ-complex? z1 z2)
    (and (equ? (complex-real-part z1) (complex-real-part z2))
         (equ? (complex-imag-part z1) (complex-imag-part z2))))
  (define (=zero-complex? x) (=zero? (complex-magnitude x)))
  (define (addd-complex z1 z2 z3)
    (make-from-real-imag (addd (complex-real-part z1)
                               (complex-real-part z2)
                               (complex-real-part z3))
                         (addd (complex-imag-part z1)
                               (complex-imag-part z2)
                               (complex-imag-part z3))))
  …
  ;; interface to rest of the system
  …
  'done)
I've included the addd procedure we added in exercise 2.82 as well for completeness.

Okay, let's see this all in action:
> (square-root (make-integer -1))
(complex rectangular (integer . 0) integer . 1)
> (make-complex-from-real-imag (make-real 3.0) (make-rational 1 2))
(complex rectangular (integer . 3) rational 1 . 2)
> (add (make-complex-from-real-imag (make-integer 4) (make-rational 2 4))
       (make-integer 4))
(complex rectangular (integer . 8) rational 1 . 2)
> (sub (make-complex-from-real-imag (make-real 1.5) (make-integer 3))
       (make-complex-from-real-imag (make-rational 1 2) (make-integer 3)))
(integer . 1)

No comments:

Post a Comment