2012-06-04

SICP Exercise 2.91: Dividing Polynomials

A univariate polynomial can be divided by another one to produce a polynomial quotient and a polynomial remainder. For example,
x5 - 1 = x3 + x, remainder x - 1
x3 - 1
Division can be performed via long division. That is, divide the highest-order term of the dividend by the highest-order term of the divisor. The result is the first term of the quotient. Next, multiply the result by the divisor, subtract that from the dividend, and produce the rest of the answer by recursively dividing the difference by the divisor. Stop when the order of the divisor exceeds the order of the dividend and declare the dividend to be the remainder. Also, if the dividend ever becomes zero, return zero as both quotient and remainder.
We can design a div-poly procedure on the model of add-poly and mul-poly. The procedure checks to see if the two polys have the same variable. If so, div-poly strips off the variable and passes the problem to div-terms, which performs the division operation on term lists. Div-poly finally reattaches the variable to the result supplied by div-terms. It is convenient to design div-terms to compute both the quotient and the remainder of a division. Div-terms can take two term lists as arguments and return a list of the quotient term list and the remainder term list.
Complete the following definition of div-terms by filling in the missing expressions. Use this to implement div-poly, which takes two polys as arguments and returns a list of the quotient and remainder polys.
(define (div-terms L1 L2)
  (if (empty-termlist? L1)
      (list (the-empty-termlist) (the-empty-termlist))
      (let ((t1 (first-term L1))
            (t2 (first-term L2)))
        (if (> (order t2) (order t1))
            (list (the-empty-termlist) L1)
            (let ((new-c (div (coeff t1) (coeff t2)))
                  (new-o (- (order t1) (order t2))))
              (let ((rest-of-result
                     <compute rest of result recursively>
                     ))
                <form complete result>
                ))))))
If you've been following along then you'll note that, due to the internal operations that are defined within of our two term-list representations (which you can see in the previous exercise), this implementation of div-terms is only applicable to sparse term-lists. We represent dense term lists as a list of the coefficients only. As a result there are no order or coeff operations, the operation first-term simply returns the coefficient of the first term rather than a term representation, and we have a different operation, term-list-order, which acts on the whole term-list and calculates the order of the highest term in the list. So for this exercise we'll start by producing a full implementation for polynomials that use sparse term-list representations only and then derive a similar implementation for dense term-lists.

Note that had we chosen a different, lower-level, interface to our term-lists then this would not be an issue. We could have had sparse and dense term-lists expose order, coeff, first-term, rest-terms, adjoin-term, the-empty-termlist and empty-termlist?. If we had done so then this would have meant that the implementations of add-terms, mul-terms, div-terms and so on would have been identical regardless of term-list representation and so could have resided in the polynomial package rather than in the two term-list packages.

However there are issues with this approach. The first-term operation would need to return a representation of a term (i.e. order and coefficient), which would lead to a lot of term creation for dense term-lists. Worse still, adjoin-term would need to take a tagged representation of a term whose tag would be stripped off when the actual operations were applied. This would mean that both the sparse and dense term-list representations would need to be able to manipulate the internal representation of a term, destroying the encapsulation. We could prevent this by having apply-generic only strip the tags off objects that are of the "type" of the operation (so, for example, only objects tagged as sparse-terms would have their tags stripped when invoking adjoin-terms for a sparse term-list).

I did consider this approach when tackling the previous exercise, as it would allow for more sharing of code. However, on examining the implications of the approach I felt it would result in overly complex changes that would lead to additional issues. Feel free to give it a go yourself though!

Anyway, on with the exercise...

The authors have provided us with a template for a recursive implementation of div-terms that includes the terminating conditions (i.e. what happens when dividend becomes zero, and what happens when the order of the divisor exceeds the order of the dividend) and the calculation of the components of the next term in the result. To complete the implementation we have to flesh out two sections of the operation: the calculation of the rest-of-result value and the production of the final result.

Calculating rest-of-results

The calculation of the rest-of-result is described above as "multiply the result by the divisor, subtract that from the dividend, and produce the rest of the answer by recursively dividing the difference by the divisor." The result here is the order (new-o) and coefficient (new-c) of the next term in the final result. To perform the multiplication we can convert the result into a term, using make-term and then invoke mul-term-by-all-terms on that term and the divisor (L2):
(mul-term-by-all-terms (make-term new-o new-c) L2)
The subtraction is slightly more tricky, as we don't actually have a subtraction operation within the sparse term-list package (nor do we have one within the dense term-list package). Such an operation is easily derived from the sub-poly implementation however: we simply negate the subtrahend (the result of the multiplication) using negate-terms, and then add this to the minuend (L1) using add-terms. This gives us the following implementation:
(add-terms L1
           (negate-terms
            (mul-term-by-all-terms (make-term new-o new-c)
                                   L2)))
The final step in this is to produce the rest-of-result by dividing the result of this calculation by the divisor (L2). This is achieved by a recursive call to div-terms, and gives us the following completed implementation of the calculation:
(div-terms 
 (add-terms L1
            (negate-terms
             (mul-term-by-all-terms (make-term new-o new-c)
                                    L2)))
 L2)

Producing the Final Result

The final result we need to produce is "a list of the quotient term list and the remainder term list," and we know that the first term of the quotient is the term formed by new-o and new-c. The rest of the quotient comes from the recursive call, and is the first item of the list that this call returns. So in order to produce the quotient we need to create the term from new-o and new-c and adjoin it to the car of rest-of-results:
(adjoin-term (make-term new-o new-c) (car rest-of-result))
The remainder has already been calculated for us - it's just the second item in the rest-of-result. So we can calculate the result as follows:
(list (adjoin-term (make-term new-o new-c)
                   (car rest-of-result))
      (cadr rest-of-result))

Putting it All Together

Okay, so we can now fill in the blanks and complete div-terms as follows:
(define (div-terms L1 L2)
  (if (empty-termlist? L1)
      (list (the-empty-termlist) (the-empty-termlist))
      (let ((t1 (first-term L1))
            (t2 (first-term L2)))
        (if (> (order t2) (order t1))
            (list (the-empty-termlist) L1)
            (let ((new-c (div (coeff t1) (coeff t2)))
                  (new-o (- (order t1) (order t2))))
              (let ((rest-of-result
                     (div-terms 
                      (add-terms L1
                                 (negate-terms
                                  (mul-term-by-all-terms (make-term new-o new-c)
                                                         L2)))
                      L2)))
                (list (adjoin-term (make-term new-o new-c)
                                   (car rest-of-result))
                      (cadr rest-of-result))))))))
You'll note that we're creating the first term in the quotient twice here: once during the calculation of rest-of-terms and once when we produce the final result. We can remove this re-creation of the term, and collapse things slightly, with the following modification:
(define (div-terms L1 L2)
  (if (empty-termlist? L1)
      (list (the-empty-termlist) (the-empty-termlist))
      (let ((t1 (first-term L1))
            (t2 (first-term L2)))
        (if (> (order t2) (order t1))
            (list (the-empty-termlist) L1)
            (let* ((new-c (div (coeff t1) (coeff t2)))
                   (new-o (- (order t1) (order t2)))
                   (new-t (make-term new-o new-c))
                   (rest-of-result
                    (div-terms
                     (add-terms L1 (negate-terms (mul-term-by-all-terms new-t L2)))
                      L2)))
                (list (adjoin-term new-t (car rest-of-result))
                      (cadr rest-of-result)))))))
Let's install this into the table of operations... As with the other arithmetic operations we need to tag the term-lists we've produced. However, we can't simply tag the result of div-terms as this is a list containing two (untagged) sparse term-lists. Instead we have to split the list apart into its two components, tag each one individually, and then put them back together. Oh, and because we're not returning a tagged type apply-generic will not automatically coerce the result to a dense term-list if that's appropriate so we'll have to do that ourselves. We can do that by using sparse-terms->dense-terms to perform the tagging:
(put 'div '(sparse-terms sparse-terms) 
     (lambda (t1 t2)
       (let ((result (div-terms t1 t2)))
         (list (sparse-terms->dense-terms (car result))
               (sparse-terms->dense-terms (cadr result))))))
This split-modify-recombine pattern repeats itself in the polynomial package...

Our div-poly operation can perform the usual work to test the variables of the two polynomials passed to it and to select the appropriate variable to add to the term-lists generated. However, similar to the tagging of the term-lists, this variable needs to be added to both term-lists generated by div-terms. As a result we need to take the result of invoking div on the numerator and denominator term-lists, split it into the quotient and remainder term-lists, add the variable to both using make-poly and then recombine them:
(define (div-poly p1 p2)
  (if (same-variable? (variable p1) (variable p2))
      (let ((variable (select-variable p1 p2))
            (result (div (term-list p1) (term-list p2))))
        (list (make-poly variable (car result))
              (make-poly variable (cadr result))))
      (error "Polys not in same var -- DIV-POLY"
             (list p1 p2))))
Finally we need to install div-poly into the table of operations. Again, we need to split, modify and recombine. The modification in this case is tagging and then dropping the result (using drop) as, again, apply-generic will not be able to do this for us automatically:
  (put 'div '(polynomial polynomial) 
       (lambda (p1 p2)
         (let ((result (div-poly p1 p2)))
           (list (drop (tag (car result)))
                 (drop (tag (cadr result)))))))

Rinse and Repeat

As we noted at the start of the exercise, due to the way in which we've implemented the sparse and dense term-list packages, the implementation provided by the authors only applies to sparse term-lists. We need to provide a separate dense term-list implementation of div-terms that uses the internal operations available to us within that package. Thankfully deriving this implementation is straightforward.

The overall form is the same, but we need to alter the implementation to account for the different representation: a dense term-list is represented as a list of coefficients rather than a list of terms. So we obtain the order of the highest term in a term-list using the term-list itself, rather than the first term, and the coefficient of the highest term is the value returned by first-term. Also we don't need to build terms - mul-term-by-all-terms and adjoin-term work directly with the order and coefficient. Here's the implementation for dense term-lists:
  (define (div-terms L1 L2)
    (if (empty-termlist? L1)
        (list (the-empty-termlist) (the-empty-termlist))
        (let ((order-1 (term-list-order L1))
              (order-2 (term-list-order L2)))
          (if (> order-2 order-1)
              (list (the-empty-termlist) L1)
              (let* ((new-c (div (first-term L1) (first-term L2)))
                     (new-o (- order-1 order-2))
                     (rest-of-result
                      (div-terms
                       (add-terms L1
                                  (negate-terms (mul-term-by-all-terms new-o new-c L2)))
                       L2)))
                (list (adjoin-term new-o new-c (car rest-of-result))
                      (cadr rest-of-result)))))))
Installing this operation in the table is similar to installing the sparse term-list operation except that instead of using the coercion procedure sparse-terms->dense-terms we use to-best-representation, which performs the coercion in the opposite direction, raising the term-list to a sparse term-list if there are enough zero terms:
  (put 'div '(dense-terms dense-terms) 
       (lambda (t1 t2)
         (let ((result (div-terms t1 t2)))
           (list (to-best-representation (car result))
                 (to-best-representation (cadr result))))))

The Proof's In the Pudding

Okay, so let's give it a spin...

First we'll define a bunch of polynomials... We'll define numerators and denominators (as both sparse and dense term-lists) that will let us perform the following calculations:
  • (x5 - 1) / (x2 - 1) = x3 + x, remainder x - 1
  • (2x2 + 2) / (x2 + 1) = 2, remainder 0
  • (3x4 + 7x3 + 6) / (0.5x4 + x3 + 3) = 6, remainder x3 - 12
Here they are defined:
(define sparse-numerator-1
  (make-polynomial-from-terms 'x
                              (list (make-term 5 (make-integer 1))
                                    (make-term 0 (make-integer -1)))))

(define sparse-denominator-1
  (make-polynomial-from-terms 'x
                              (list (make-term 2 (make-integer 1))
                                    (make-term 0 (make-integer -1)))))

(define sparse-numerator-2
  (make-polynomial-from-terms 'x
                              (list (make-term 2 (make-integer 2))
                                    (make-term 0 (make-integer 2)))))

(define sparse-denominator-2
  (make-polynomial-from-terms 'x
                              (list (make-term 2 (make-integer 1))
                                    (make-term 0 (make-integer 1)))))

(define sparse-numerator-3
  (make-polynomial-from-terms 'x
                              (list (make-term 4 (make-integer 3))
                                    (make-term 3 (make-integer 7))
                                    (make-term 0 (make-integer 6)))))

(define sparse-denominator-3
  (make-polynomial-from-terms 'x
                              (list (make-term 4 (make-real 0.5))
                                    (make-term 3 (make-integer 1))
                                    (make-term 0 (make-integer 3)))))

(define dense-numerator-1
  (make-polynomial-from-coeffs 'x
                               (list (make-integer 1)
                                     zero
                                     zero
                                     zero
                                     zero
                                     (make-integer -1))))

(define dense-denominator-1
  (make-polynomial-from-coeffs 'x
                               (list (make-integer 1)
                                     zero
                                     (make-integer -1))))

(define dense-numerator-2
  (make-polynomial-from-coeffs 'x
                               (list (make-integer 2)
                                     zero
                                     (make-integer 2))))

(define dense-denominator-2
  (make-polynomial-from-coeffs 'x
                               (list (make-integer 1)
                                     zero
                                     (make-integer 1))))

(define dense-numerator-3
  (make-polynomial-from-coeffs 'x
                               (list (make-integer 3)
                                     (make-integer 7)
                                     zero
                                     zero
                                     (make-integer 6))))

(define dense-denominator-3
  (make-polynomial-from-coeffs 'x
                               (list (make-real 0.5)
                                     (make-integer 1)
                                     zero
                                     zero
                                     (make-integer 3))))
And here they are evaluated:
> (div sparse-numerator-1 sparse-denominator-1)
((polynomial x sparse-terms (term 3 (integer . 1)) (term 1 (integer . 1)))
 (polynomial x dense-terms (integer . 1) (integer . -1)))
> (div sparse-numerator-2 sparse-denominator-2)
((integer . 2)
 (integer . 0))
> (div sparse-numerator-3 sparse-denominator-3)
((integer . 6)
 (polynomial x sparse-terms (term 3 (integer . 1)) (term 0 (integer . -12))))
> (div dense-numerator-1 dense-denominator-1)
((polynomial x sparse-terms (term 3 (integer . 1)) (term 1 (integer . 1)))
 (polynomial x dense-terms (integer . 1) (integer . -1)))
> (div dense-numerator-2 dense-denominator-2)
((integer . 2)
 (integer . 0))
> (div dense-numerator-3 dense-denominator-3)
((integer . 6)
 (polynomial x sparse-terms (term 3 (integer . 1)) (term 0 (integer . -12))))
> (div dense-numerator-1 sparse-denominator-1)
((polynomial x sparse-terms (term 3 (integer . 1)) (term 1 (integer . 1)))
 (polynomial x dense-terms (integer . 1) (integer . -1)))
> (div sparse-numerator-1 dense-denominator-1)
((polynomial x sparse-terms (term 3 (integer . 1)) (term 1 (integer . 1)))
 (polynomial x dense-terms (integer . 1) (integer . -1)))

2012-05-31

SICP Exercise 2.90: Supporting Dense and Sparse Polynomials - Part 3

Suppose we want to have a polynomial system that is efficient for both sparse and dense polynomials. One way to do this is to allow both kinds of term-list representations in our system. The situation is analogous to the complex-number example of section 2.4, where we allowed both rectangular and polar representations. To do this we must distinguish different types of term lists and make the operations on term lists generic. Redesign the polynomial system to implement this generalization. This is a major effort, not a local change.
Okay, we've talked about it at length and updated the system so that polynomials form a part of our overall generic arithmetic system. Let's put it to the test and check it all works.

First we should define a few polynomials to work with:
(define dense
        (make-polynomial-from-coeffs 'x
                                     (list (make-integer 4)
                                           (make-integer 3)
                                           (make-integer 2)
                                           (make-integer 1)
                                           zero)))

(define dense-with-many-zeros
        (make-polynomial-from-coeffs 'x
                                     (list (make-integer 42)
                                           zero
                                           zero
                                           zero
                                           zero
                                           zero
                                           (make-integer -1))))

(define sparse
        (make-polynomial-from-terms 'x 
                                    (list (make-term 5 (make-integer 5))
                                          (make-term 3 (make-integer 3))
                                          (make-term 1 (make-integer 1)))))

(define another-sparse
        (make-polynomial-from-terms 'x
                                    (list (make-term 5 (make-integer 5))
                                          (make-term 3 (make-integer 3))
                                          (make-term 1 (make-integer 1))
                                          (make-term 0 (make-integer 3)))))

(define very-sparse
        (make-polynomial-from-terms 'x
                                    (list (make-term 50 (make-integer 150))
                                          (make-term 10 (make-integer 11))
                                          (make-term 0 (make-integer 1)))))

(define polypoly
        (make-polynomial-from-coeffs
         'x
         (list (make-polynomial-from-coeffs 'y
                                            (list (make-integer 2)
                                                  (make-integer 1))))))
Now let's try out a bunch of operations and see what happens:
> (add polypoly dense)
(polynomial x
            dense-terms
            (integer . 4)
            (integer . 3)
            (integer . 2)
            (integer . 1)
            (polynomial y
                        dense-terms
                        (integer . 2)
                        (integer . 1)))
> (add polypoly polypoly)
(polynomial x
            dense-terms
            (polynomial y
                        dense-terms
                        (integer . 4)
                        (integer . 2)))
> (add (add polypoly polypoly) (make-integer 3))
(polynomial x
            dense-terms
            (polynomial y
                        dense-terms
                        (integer . 4)
                        (integer . 5)))
> (add dense dense-with-many-zeros)
(polynomial x
            dense-terms
            (integer . 42)
            (integer . 0)
            (integer . 4)
            (integer . 3)
            (integer . 2)
            (integer . 1)
            (integer . -1))
> (add dense-with-many-zeros dense-with-many-zeros)
(polynomial x
            sparse-terms
            (term 6 (integer . 84))
            (term 0 (integer . -2)))
> (add sparse sparse)
(polynomial x
            sparse-terms
            (term 5 (integer . 10))
            (term 3 (integer . 6))
            (term 1 (integer . 2)))
> (add sparse another-sparse)
(polynomial x
            sparse-terms
            (term 5 (integer . 10))
            (term 3 (integer . 6))
            (term 1 (integer . 2))
            (term 0 (integer . 3)))
> (add very-sparse sparse)
(polynomial x
            sparse-terms
            (term 50 (integer . 150))
            (term 10 (integer . 11))
            (term 5 (integer . 5))
            (term 3 (integer . 3))
            (term 1 (integer . 1))
            (term 0 (integer . 1)))
> (mul sparse dense)
(polynomial x
            sparse-terms
            (term 9 (integer . 20))
            (term 8 (integer . 15))
            (term 7 (integer . 22))
            (term 6 (integer . 14))
            (term 5 (integer . 10))
            (term 4 (integer . 6))
            (term 3 (integer . 2))
            (term 2 (integer . 1)))
> (add dense sparse)
(polynomial x
            dense-terms
            (integer . 5)
            (integer . 4)
            (integer . 6)
            (integer . 2)
            (integer . 2)
            (integer . 0))
> (sub sparse dense)
(polynomial x
            sparse-terms
            (term 5 (integer . 5))
            (term 4 (integer . -4))
            (term 2 (integer . -2)))
> (negate very-sparse)
(polynomial x
            sparse-terms
            (term 50 (integer . -150))
            (term 10 (integer . -11))
            (term 0 (integer . -1)))
> (sub (add dense (make-integer 1)) dense)
(integer . 1)

SICP Exercise 2.90: Supporting Dense and Sparse Polynomials - Part 2

Suppose we want to have a polynomial system that is efficient for both sparse and dense polynomials. One way to do this is to allow both kinds of term-list representations in our system. The situation is analogous to the complex-number example of section 2.4, where we allowed both rectangular and polar representations. To do this we must distinguish different types of term lists and make the operations on term lists generic. Redesign the polynomial system to implement this generalization. This is a major effort, not a local change.
We've already discussed how we're going to approach this, so let's get onto the code. For a change I'm going to post all of the code involved in the system so far, rather than just the changes, so I'll include all the code for the integer, rational, real and complex packages too, along with all the supporting framework. We've been working on the system for quite a while now, so this'll help to provide a check point for where we've got to.

I'll intersperse the code with commentary where I feel it's appropriate... Oh, and it's worth noting that I've tweaked and tidied code from previous exercises as I've gone along, so some bits won't correspond exactly to what's gone before.

Here we go...
; If you're using DrRacket you might need to uncomment these...
;(require rnrs/base-6)
;(require rnrs/mutable-pairs-6)


;;;
;;; Basic definitions
;;;

(define false #f)
(define true #t)


;;;
;;; Look-up table from Creating local tables in section 3.3.3
;;;

(define (assoc key records)
  (cond ((null? records) false)
        ((equal? key (caar records)) (car records))
        (else (assoc key (cdr records)))))

(define (make-table)
  (let ((local-table (list '*table*)))
    (define (lookup key-1 key-2)
      (let ((subtable (assoc key-1 (cdr local-table))))
        (if subtable
            (let ((record (assoc key-2 (cdr subtable))))
              (if record
                  (cdr record)
                  false))
            false)))
    (define dump local-table)
    (define (insert! key-1 key-2 value)
      (let ((subtable (assoc key-1 (cdr local-table))))
        (if subtable
            (let ((record (assoc key-2 (cdr subtable))))
              (if record
                  (set-cdr! record value)
                  (set-cdr! subtable
                            (cons (cons key-2 value)
                                  (cdr subtable)))))
            (set-cdr! local-table
                      (cons (list key-1
                                  (cons key-2 value))
                            (cdr local-table)))))
      'ok)    
    (define (dispatch m)
      (cond ((eq? m 'lookup-proc) lookup)
            ((eq? m 'insert-proc!) insert!)
            ((eq? m 'dump) dump)
            (else (error "Unknown operation -- TABLE" m))))
    dispatch))

(define operation-table (make-table))
(define get (operation-table 'lookup-proc))
(define put (operation-table 'insert-proc!))
(define dump-table (operation-table 'dump))

(define (put-coercion source-type target-type proc)
  (put 'coercion (list source-type target-type) proc))

(define (get-coercion source-type target-type)
  (get 'coercion (list source-type target-type)))


;;;
;;; Tagging utils
;;;

(define (attach-tag type-tag contents)
  (cons type-tag contents))
(define (type-tag datum)
  (if (pair? datum)
      (car datum)
      (error "Bad tagged datum -- TYPE-TAG" datum)))
(define (contents datum)
  (if (pair? datum)
      (cdr datum)
      (error "Bad tagged datum -- CONTENTS" datum)))


;;;
;;; Generic Arithmetic operations
;;;

(define (add x y) (apply-generic 'add x y))
(define (sub x y) (apply-generic 'sub x y))
(define (mul x y) (apply-generic 'mul x y))
(define (div x y) (apply-generic 'div x y))
(define (equ? x y) (apply-generic 'equ? x y))
(define (=zero? x) (apply-generic '=zero? x))
(define (addd x y z) (apply-generic 'addd x y z))
(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))
(define (negate x) (apply-generic 'negate x))
First item of note: we update the tower-of-types to include the polynomial package. We also add in the dense-terms and sparse-terms packages to allow coercion between these two types, but we'll not install coercion procedures for polynomial to dense-terms or vice versa, effectively splitting the tower-of-types into two towers.
;;;
;;; Tower-of-types
;;;

(define tower-of-types '(integer rational real complex polynomial
                         dense-terms sparse-terms))


;;;
;;; Generic Operation Framework
;;;

(define (find-highest-type list-of-types)
  (define (filter-type type list-to-filter)
    (cond ((null? list-to-filter) '())
          ((eq? (car list-to-filter) type) (filter-type type (cdr list-to-filter)))
          (else (cons (car list-to-filter) (filter-type type (cdr list-to-filter))))))
  (define (find-highest highest remaining-tower remaining-list)
    (cond ((null? remaining-list) highest)
          ((null? remaining-tower)
           (error "Cannot find highest type from non-tower types -- FIND-HIGHEST-TYPE"
                  remaining-list))
          (else (find-highest (car remaining-tower)
                              (cdr remaining-tower)
                              (filter-type (car remaining-tower) remaining-list)))))
  (find-highest #f tower-of-types list-of-types))
  
(define (raise-to type value)
  (cond ((eq? type (type-tag value)) value)
        ((memq type tower-of-types) (raise-to type (raise value)))
        (else (error "Cannot raise to non-tower type -- RAISE-TO"
                     (list type tower-of-types)))))

(define (raise-all-to type values)
  (if (null? values)
      '()
      (cons (raise-to type (car values)) (raise-all-to type (cdr values)))))

As noted in the preceding discussion, in order to allow a missing coercion operation to indicate a point at which a type should not be lowered any further (and so to allow multiple towers-of-types to be represented by a single tower) we modify apply-raise so that it returns the value passed to it, rather than raising an error, if it can't find an appropriate coercion operation.
(define (apply-raise x types)
  (cond ((null? types)
         (error "Type not found in the tower-of-types"
                (list (type-tag x) tower-of-types)))
        ((eq? (type-tag x) (car types))
         (if (null? (cdr types))
             x
             (let ((raiser (get-coercion (type-tag x) (cadr types))))
               (if raiser
                   (raiser (contents x))
                   x))))
        (else (apply-raise x (cdr types)))))

(define (raise x)
  (apply-raise x tower-of-types))

(define (project x)
  (apply-raise x (reverse tower-of-types)))

(define (drop x)
  (let ((dropped (project x)))
    (if (eq? (type-tag x) (type-tag dropped))
        x
        (let ((raised (raise dropped)))
          (if (not (equ? x raised))
              x
              (drop dropped))))))

(define (raise-to type value)
  (cond ((eq? type (type-tag value)) value)
        ((memq type tower-of-types) (raise-to type (raise value)))
        (else (error "Cannot raise to non-tower type -- RAISE-TO"
                     (list type tower-of-types)))))

(define (in-tower? value)
  (and (pair? value) (memq (type-tag value) tower-of-types)))

(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)))))

(define (apply-generic op . args)
  (define (find-and-apply-op)
    (let* ((type-tags (map type-tag args))
           (proc (get op type-tags)))
      (if proc
          (apply proc (map contents args))
          (if (> (length args) 1)
              (let* ((highest-type (find-highest-type type-tags))
                     (mapped-args (raise-all-to highest-type args))
                     (mapped-types (map type-tag mapped-args))
                     (mapped-proc (get op mapped-types)))
                (if mapped-proc
                    (apply mapped-proc (map contents mapped-args))
                    (error
                     "No method for these types -- APPLY-GENERIC"
                     (list op type-tags))))))))
  (let ((result (find-and-apply-op)))
    (if (in-tower? result)
        (drop result)
        result)))
The standard numerical type representations follow. They're mostly unchanged apart from the addition of support for the negate operation, and the addition of a coercion operation to raise a complex type to a polynomial.
;;;
;;; Integer package
;;;

(define (install-integer-package)
  (define (tag x)
    (attach-tag 'integer x))    
  (define (integer->rational i) (make-rational i 1))
  (put 'add '(integer integer)
       (lambda (x y) (tag (+ x y))))
  (put 'sub '(integer integer)
       (lambda (x y) (tag (- x y))))
  (put 'mul '(integer integer)
       (lambda (x y) (tag (* x y))))
  (put 'div '(integer integer)
       (lambda (x y) (make-rational x y)))
  (put 'equ? '(integer integer) =)
  (put '=zero? '(integer)
       (lambda (x) (= 0 x)))
  (put 'addd '(integer integer integer)
       (lambda (x y z) (tag (+ x y z))))
  (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))))
  (put 'negate '(integer)
       (lambda (x) (tag (- x))))
  (put 'make 'integer
       (lambda (x) (if (integer? x)
                       (tag x)
                       (error "non-integer value" x))))
  (put-coercion 'integer 'rational integer->rational)
  'done)

(define (make-integer n)
  ((get 'make 'integer) n))


;;;
;;; Rational numbers package
;;;

(define (install-rational-package)
  ;; internal procedures
  (define (numer x) (car x))
  (define (denom x) (cdr x))
  (define (make-rat n d)
    (if (and (integer? n) (integer? d))
        (let ((g (gcd n d)))
          (cons (/ n g) (/ d g)))
        (error "non-integer numerator or denominator"
               (list n d))))
  (define (add-rat x y)
    (make-rat (+ (* (numer x) (denom y))
                 (* (numer y) (denom x)))
              (* (denom x) (denom y))))
  (define (sub-rat x y)
    (make-rat (- (* (numer x) (denom y))
                 (* (numer y) (denom x)))
              (* (denom x) (denom y))))
  (define (mul-rat x y)
    (make-rat (* (numer x) (numer y))
              (* (denom x) (denom y))))
  (define (div-rat x y)
    (make-rat (* (numer x) (denom y))
              (* (denom x) (numer y))))
  (define (equ? x y)
    (and (= (numer x) (numer y))
         (= (denom x) (denom y))))
  (define (=zero? x) (= (numer x) 0))
  (define (addd x y z)
    (make-rat (+ (* (numer x) (denom y) (denom z))
                 (* (denom x) (numer y) (denom z))
                 (* (denom x) (denom y) (numer z)))
              (* (denom x) (denom y) (denom z))))
  (define (rational->real r) (make-real (/ (numer r) (denom r))))
  (define (rational->integer r) (make-integer (round (/ (numer r) (denom r)))))
  (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))))
  (define (negate-rat x)
    (make-rat (- (numer x)) (denom x)))
  ;; interface to rest of the system
  (define (tag x) (attach-tag 'rational x))
  (put 'add '(rational rational)
       (lambda (x y) (tag (add-rat x y))))
  (put 'sub '(rational rational)
       (lambda (x y) (tag (sub-rat x y))))
  (put 'mul '(rational rational)
       (lambda (x y) (tag (mul-rat x y))))
  (put 'div '(rational rational)
       (lambda (x y) (tag (div-rat x y))))
  (put 'equ? '(rational rational) equ?)
  (put '=zero? '(rational) =zero?)
  (put 'addd '(rational rational rational)
       (lambda (x y z) (tag (addd x y z))))
  (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))))
  (put 'negate '(rational)
       (lambda (x) (tag (negate-rat x))))
  (put 'make 'rational
       (lambda (n d) (tag (make-rat n d))))
  (put-coercion 'rational 'real rational->real)
  (put-coercion 'rational 'integer rational->integer)
  'done)

(define (make-rational n d)
  ((get 'make 'rational) n d))


;;;
;;; Real package
;;;

(define (install-real-package)
  (define (tag x)
    (attach-tag 'real x))    
  (define (real->complex r) (make-complex-from-real-imag (tag r) (tag 0)))
  (define (real->rational r) (make-rational (inexact->exact (numerator r))
                                            (inexact->exact (denominator r))))
  (put 'add '(real real)
       (lambda (x y) (tag (+ x y))))
  (put 'sub '(real real)
       (lambda (x y) (tag (- x y))))
  (put 'mul '(real real)
       (lambda (x y) (tag (* x y))))
  (put 'div '(real real)
       (lambda (x y) (tag (/ x y))))
  (put 'equ? '(real real) =)
  (put '=zero? '(real)
       (lambda (x) (= 0 x)))
  (put 'addd '(real real real)
       (lambda (x y z) (tag (+ x y z))))
  (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))))
  (put 'negate '(real)
       (lambda (x) (tag (- x))))
  (put 'make 'real
       (lambda (x) (if (real? x)
                       (tag x)
                       (error "non-real value" x))))
  (put-coercion 'real 'complex real->complex)
  (put-coercion 'real 'rational real->rational)
  'done)

(define (make-real n)
  ((get 'make 'real) n))


;;;
;;; Complex (rectangular) number package
;;;

(define (install-rectangular-package)
  ;; internal procedures
  (define (real-part z) (car z))
  (define (imag-part z) (cdr z))
  (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
  (define (tag x) (attach-tag 'rectangular x))
  (put 'real-part '(rectangular) real-part)
  (put 'imag-part '(rectangular) imag-part)
  (put 'magnitude '(rectangular) magnitude)
  (put 'angle '(rectangular) angle)
  (put 'make-from-real-imag 'rectangular 
       (lambda (x y) (tag (make-from-real-imag x y))))
  (put 'make-from-mag-ang 'rectangular 
       (lambda (r a) (tag (make-from-mag-ang r a))))
  'done)

;;;
;;; Complex (polar) number package
;;;

(define (install-polar-package)
  ;; internal procedures
  (define (magnitude z) (car z))
  (define (angle z) (cdr z))
  (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
  (define (tag x) (attach-tag 'polar x))
  (put 'real-part '(polar) real-part)
  (put 'imag-part '(polar) imag-part)
  (put 'magnitude '(polar) magnitude)
  (put 'angle '(polar) angle)
  (put 'make-from-real-imag 'polar
       (lambda (x y) (tag (make-from-real-imag x y))))
  (put 'make-from-mag-ang 'polar 
       (lambda (r a) (tag (make-from-mag-ang r a))))
  'done)


;;;
;;; Generic complex number procedures
;;;

;;; N.B. Name-clash with complex number support in rnrs/base-6
;;; - prefixing complex ops with complex-

(define (complex-real-part z) (apply-generic 'real-part z))
(define (complex-imag-part z) (apply-generic 'imag-part z))
(define (complex-magnitude z) (apply-generic 'magnitude z))
(define (complex-angle z) (apply-generic 'angle z))
Note the body of the complex->polynomial coercion operation below: (make-zero-order-polynomial-from-coeff (drop (tag z))). We re-tag the complex number as such, as it will have no tag at this point, and then drop it as far as possible before using a special creation operation. The former is necessary in order to simplify the type as much as possible before creating the polynomial (and so to undo any coercion that may have occurred getting to this level in the tree). The latter allows us to encapsulate the special 'unbound marker within the polynomial package.
;;;
;;; Complex number package
;;;

(define (install-complex-package)
  ;; imported procedures from rectangular and polar packages
  (define (make-from-real-imag x y)
    ((get 'make-from-real-imag 'rectangular) x y))
  (define (make-from-mag-ang r a)
    ((get 'make-from-mag-ang 'polar) r a))
  ;; 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))))
  (define (complex->real z) (raise-to 'real (complex-real-part z)))
  (define (complex->polynomial z) (make-zero-order-polynomial-from-coeff (drop (tag z))))
  (define (negate-complex z)
    (make-from-real-imag (negate (complex-real-part z))
                         (negate (complex-imag-part z))))
  ;; interface to rest of the system
  (define (tag z) (attach-tag 'complex z))
  (put 'add '(complex complex)
       (lambda (z1 z2) (tag (add-complex z1 z2))))
  (put 'sub '(complex complex)
       (lambda (z1 z2) (tag (sub-complex z1 z2))))
  (put 'mul '(complex complex)
       (lambda (z1 z2) (tag (mul-complex z1 z2))))
  (put 'div '(complex complex)
       (lambda (z1 z2) (tag (div-complex z1 z2))))
  (put 'equ? '(complex complex) equ-complex?)
  (put '=zero? '(complex) =zero-complex?)
  (put 'addd '(complex complex complex)
       (lambda (z1 z2 z3) (tag (addd-complex z1 z2 z3))))
  ;;; N.B. Name-clash with complex number support in rnrs/base-6
  ;;; - prefixing complex ops with complex-
  (put 'real-part '(complex) complex-real-part)
  (put 'imag-part '(complex) complex-imag-part)
  (put 'magnitude '(complex) complex-magnitude)
  (put 'angle '(complex) complex-angle)
  (put 'negate '(complex)
       (lambda (z) (tag (negate-complex z))))
  (put 'make-from-real-imag 'complex
       (lambda (x y) (tag (make-from-real-imag x y))))
  (put 'make-from-mag-ang 'complex
       (lambda (r a) (tag (make-from-mag-ang r a))))
  (put-coercion 'complex 'real complex->real)
  (put-coercion 'complex 'polynomial complex->polynomial)
  'done)

(define (make-complex-from-real-imag x y)
  ((get 'make-from-real-imag 'complex) x y))
(define (make-complex-from-mag-ang r a)
  ((get 'make-from-mag-ang 'complex) r a))


;;;
;;; Install numerical type packages
;;;

(install-rectangular-package)
(install-polar-package)

(install-integer-package)
(install-rational-package)
(install-real-package)
(install-complex-package)
Zero is a useful concept... We used it in exercise 2.89 to avoid constantly recreating zero coefficients. We'll use it again here!
;;;
;;; Zero
;;;
(define zero (make-integer 0))
We're going to be able to construct both sparse and dense term-lists using terms, so we pull out the concept of a term into its own package. While it uses generic operations to encapsulate its internal representation, it doesn't exist within the tower-of-types.
;;;
;;; Terms
;;;

(define (install-term-package)
  ;; internal procedures
  ;; representation of terms
  (define (make-term order coeff) (list order coeff))
  (define (order-term term) (car term))
  (define (coeff-term term) (cadr term))
  (define (equ-term? t1 t2)
    (and (= (order-term t1) (order-term t2))
         (equ? (coeff-term t1) (coeff-term t2))))

  ;; interface to rest of the system
  (define (tag t) (attach-tag 'term t))
  (put 'order '(term) order-term)
  (put 'coeff '(term) coeff-term)
  (put 'equ? '(term term) equ-term?)
  (put 'make 'term
       (lambda (order coeff) (tag (make-term order coeff))))
  'done)

(define (make-term order coeff)
  ((get 'make 'term) order coeff))

(define (order term)
  (apply-generic 'order term))

(define (coeff term)
  (apply-generic 'coeff term))
Before we get into the specific implementations of the term-lists, we need to know when to use which representation. As there are only two available representations in our system at the moment this is a binary choice, so we produce a predicate that can tell us when to use the sparse term-list representation, and then fall back to using the dense term-list representation in all other cases. Here's the predicate:
;;;
;;; Sparse/dense representation selector
;;;
(define (store-as-sparse? highest-order zero-terms)
  (if (>= highest-order 10)
      (> (/ zero-terms highest-order) 0.1)
      (> zero-terms (/ highest-order 5))))
The sparse term-list package is mostly derived from the implementation presented in section 2.5.3 and our work in exercises 2.87 and 2.88, with operations related to polynomial manipulation stripped out. There are some slight tweaks to the behaviour, of course...
  • We ensure that we always have at least one term in any polynomial, even if its a single order 0 term with a coefficient of 0. This makes some manipulations easier, especially coercing from polynomial to complex, as we don't have to check especially for an empty term list. We use this same trick with the dense term-list package too.
  • We add a creation operation that will create sparse term-lists from coefficients. This operation is mostly implemented by convert-to-term-list, which assigns an order to each coefficient based upon its position in the list.
  • There are various additions in order to allow coercing a sparse term-list to a dense term-list. These include operations to support testing whether or not it's valid to represent a sparse term-list as a dense term-list as well as the coercion procedure itself, which simply calls through to the appropriate dense term-list creation operation if it's a viable coercion. An implementation of the equ? operation is also added, as drop coercions need to be able to test whether projecting and then re-raising a sparse polynomial gives a value equal to the original before they can allow drops to succeed.
  • There are also various operations added to support the coercion of polynomial types to complex. The implementation of equ? is also used here, and a new operation, get-constant, is added, which returns the coefficient from the order 0 term of the term-list (or returns zero if there is no such term in the representation).
  • The interface to the rest of the system changes, so that it exposes operations on term-lists rather than polynomials, and we also expose a few more operations.
;;;
;;; Sparse term-lists
;;; 
(define (install-sparse-terms-package)
  ;; internal procedures
  ;; representation of terms and term lists
  (define (adjoin-term term term-list)
    (if (=zero? (coeff term))
        term-list
        (cons term term-list)))
  (define (the-empty-termlist) '())
  (define (first-term term-list) (car term-list))
  (define (rest-terms term-list) (cdr term-list))
  (define (empty-termlist? term-list) (null? term-list))

  ;; Ensures a term list always has at least one term
  (define (ensure-valid-term-list terms)
    (if (empty-termlist? terms)
        (list (make-term 0 zero))
        terms))
               
  ;; creation
  (define (insert-term term terms)
    (if (empty-termlist? terms)
        (adjoin-term term (the-empty-termlist))
        (let* ((head (first-term terms))
               (head-order (order head))
               (term-order (order term)))
          (cond ((> term-order head-order) (adjoin-term term terms))
                ((= term-order head-order)
                 (adjoin-term (make-term term-order (add (coeff term) (coeff head)))
                              (rest-terms terms)))
                (else (adjoin-term head (insert-term term (rest-terms terms))))))))
  (define (build-terms terms result)
    (if (null? terms)
        result
        (build-terms (cdr terms) (insert-term (car terms) result))))
  (define (make-from-terms terms)
    (build-terms terms (the-empty-termlist)))

  (define (convert-to-term-list coeffs)
    (if (null? coeffs)
        (the-empty-termlist)
        (adjoin-term (make-term (- (length coeffs) 1) (car coeffs))
                     (convert-to-term-list (cdr coeffs)))))
  (define (make-from-coeffs coeffs)
    (convert-to-term-list coeffs))

  ;; procedures used by add-poly
  (define (add-terms L1 L2)
    (cond ((empty-termlist? L1) L2)
          ((empty-termlist? L2) L1)
          (else
           (let ((t1 (first-term L1)) (t2 (first-term L2)))
             (cond ((> (order t1) (order t2))
                    (adjoin-term
                     t1 (add-terms (rest-terms L1) L2)))
                   ((< (order t1) (order t2))
                    (adjoin-term
                     t2 (add-terms L1 (rest-terms L2))))
                   (else
                    (adjoin-term
                     (make-term (order t1)
                                (add (coeff t1) (coeff t2)))
                     (add-terms (rest-terms L1)
                                (rest-terms L2)))))))))

  ;; procedures used by mul-poly
  (define (mul-terms L1 L2)
    (if (empty-termlist? L1)
        (the-empty-termlist)
        (add-terms (mul-term-by-all-terms (first-term L1) L2)
                   (mul-terms (rest-terms L1) L2))))
  (define (mul-term-by-all-terms t1 L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let ((t2 (first-term L)))
          (adjoin-term
           (make-term (+ (order t1) (order t2))
                      (mul (coeff t1) (coeff t2)))
           (mul-term-by-all-terms t1 (rest-terms L))))))

  ;; procedures used by div-poly
  (define (div-terms L1 L2)
    (if (empty-termlist? L1)
        (list (the-empty-termlist) (the-empty-termlist))
        (let ((t1 (first-term L1))
              (t2 (first-term L2)))
          (if (> (order t2) (order t1))
              (list (the-empty-termlist) L1)
              (let* ((new-c (div (coeff t1) (coeff t2)))
                     (new-o (- (order t1) (order t2)))
                     (new-t (make-term new-o new-c))
                     (rest-of-result
                      (div-terms
                       (add-terms L1
                                  (negate-terms (mul-term-by-all-terms new-t L2)))
                        L2)))
                  (list (adjoin-term new-t (car rest-of-result))
                        (cadr rest-of-result)))))))

  ;; zero test
  (define (=zero-all-terms? L)
    (cond ((empty-termlist? L) #t)
          ((not (=zero? (coeff (first-term L)))) #f)
          (else (=zero-all-terms? (rest-terms L)))))
        
  ;; Negation
  (define (negate-terms L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let ((term (first-term L)))
          (adjoin-term (make-term (order term)
                                  (negate (coeff term)))
                       (negate-terms (rest-terms L))))))

  ;; Equality
  (define (equ-terms? L1 L2)
    (cond ((empty-termlist? L1) (empty-termlist? L2))
          ((empty-termlist? L2) #f)
          (else (and (equ? (first-term L1) (first-term L2))
                     (equ-terms? (rest-terms L1) (rest-terms L2))))))

  ;; Constant
  (define (get-sparse-constant L)
    (cond ((empty-termlist? L) zero)
          ((= (order (first-term L)) 0) (coeff (first-term L)))
          (else (get-sparse-constant (rest-terms L)))))
  
  ;; Coercion
  (define (calculate-zero-terms first rest)
    (if (empty-termlist? rest)
        (order first)
        (let ((next (first-term rest)))
          (+ (- (order first) (order next) 1)
             (calculate-zero-terms next (rest-terms rest))))))
         
  (define (keep-as-sparse? L)
    (if (empty-termlist? L)
        #f
        (let ((highest-order (order (first-term L)))
              (zero-terms (calculate-zero-terms (first-term L) (rest-terms L))))
          (store-as-sparse? highest-order zero-terms))))

  (define (sparse-terms->dense-terms L)
    (if (keep-as-sparse? L)
        (tag L)
        ((get 'make-from-terms 'dense-terms) L)))

  ;; interface to rest of the system
  (define (tag t) (attach-tag 'sparse-terms (ensure-valid-term-list t)))
  (put 'add '(sparse-terms sparse-terms) 
       (lambda (t1 t2) (tag (add-terms t1 t2))))
  (put 'mul '(sparse-terms sparse-terms) 
       (lambda (t1 t2) (tag (mul-terms t1 t2))))
  (put 'div '(sparse-terms sparse-terms) 
       (lambda (t1 t2)
         (let ((result (div-terms t1 t2)))
           (list (sparse-terms->dense-terms (car result))
                 (sparse-terms->dense-terms (cadr result))))))
  (put 'equ? '(sparse-terms sparse-terms) equ-terms?)
  (put '=zero? '(sparse-terms) =zero-all-terms?)
  (put 'negate '(sparse-terms)
       (lambda (t) (tag (negate-terms t))))
  (put 'get-constant '(sparse-terms) get-sparse-constant)
  (put 'make-from-terms 'sparse-terms
       (lambda (terms) (tag (make-from-terms terms))))
  (put 'make-from-coeffs 'sparse-terms
       (lambda (coeffs) (tag (make-from-coeffs coeffs))))
  (put-coercion 'sparse-terms 'dense-terms sparse-terms->dense-terms)
  'done)
Similarly, the dense term-list package is mostly derived from the our work in exercise 2.89, with operations related to polynomial manipulation stripped out and a similar bunch of slight tweaks to the behaviour...
  • We ensure that we always have at least one term in any polynomial.
  • We add a creation operation that will create dense term-lists from terms. This operation is derived from the sparse term-list's equivalent creation operation.
  • Coercing a dense term-list to a sparse term-list is much simpler than the reverse - we add a coercion operation which simply calls through to the appropriate sparse term-list creation. However, we won't get automatic raising to sparse term-lists when a dense term-list representation is inappropriate, so we add an operation, to-best-representation, which takes the result of an arithmetic operation and either tags it as a dense term-list or raises it to a sparse term-list as appropriate. We then modify all of the arithmetic operations to use this instead of tag.
  • The same set of operations are added to support the coercion of polynomial types to complex as were added for sparse term-lists. I.e. implementations of equ? and get-constant.
  • The interface to the rest of the system changes in an identical manner to sparse term-lists.
;;;
;;; Dense term-lists
;;; 
(define (install-dense-terms-package)
  ;; internal procedures
  ;; representation of term lists
  (define (adjoin-term term-order term-coeff term-list)
    (cond ((=zero? term-coeff) term-list)
          ((= term-order (+ 1 (term-list-order term-list)))
           (cons term-coeff term-list))
          ((> term-order (term-list-order term-list))
           (adjoin-term term-order
                        term-coeff
                        (cons zero term-list)))
          (else (error "Cannot adjoin term of lower order than term list -- ADJOIN-TERM"
                       (list term-order term-coeff term-list)))))
  (define (the-empty-termlist) '())
  (define (first-term term-list) (car term-list))
  (define (rest-terms term-list) (cdr term-list))
  (define (empty-termlist? term-list) (null? term-list))
  (define (term-list-order term-list)
    (- (length term-list) 1))

  ;; Ensures a term list always has at least one term
  (define (ensure-valid-term-list terms)
    (if (empty-termlist? terms)
        (list zero)
        terms))

  ;; Creation
  (define (strip-leading-zeros coeffs)
    (cond ((empty-termlist? coeffs) (the-empty-termlist))
          ((not (=zero? (first-term coeffs))) coeffs)
          (else (make-from-coeffs (rest-terms coeffs)))))
  (define (make-from-coeffs coeffs) coeffs)
  
  (define (insert-term term terms)
    (if (empty-termlist? terms)
        (adjoin-term (order term) (coeff term) (the-empty-termlist))
        (let ((head-order (term-list-order terms))
              (term-order (order term)))
          (cond ((> term-order head-order)
                 (adjoin-term term-order (coeff term) terms))
                ((= term-order head-order)
                 (adjoin-term term-order
                              (add (coeff term) (car terms))
                              (rest-terms terms)))
                (else (adjoin-term head-order
                                   (car terms)
                                   (insert-term term (rest-terms terms))))))))
  (define (build-terms terms result)
    (if (null? terms)
        result
        (build-terms (cdr terms) (insert-term (car terms) result))))
  (define (make-from-terms terms)
    (build-terms terms (the-empty-termlist)))
  
  ;; procedures used by add-poly
  (define (add-terms L1 L2)
    (cond ((empty-termlist? L1) L2)
          ((empty-termlist? L2) L1)
          (else
           (let* ((t1 (first-term L1))
                  (t2 (first-term L2))
                  (o1 (term-list-order L1))
                  (o2 (term-list-order L2)))
             (cond ((> o1 o2)
                    (adjoin-term
                     o1 t1 (add-terms (rest-terms L1) L2)))
                   ((< o1 o2)
                    (adjoin-term
                     o2 t2 (add-terms L1 (rest-terms L2))))
                   (else
                    (adjoin-term
                     o1
                     (add t1 t2)
                     (add-terms (rest-terms L1)
                                (rest-terms L2)))))))))

  ;; procedures used by mul-poly
  (define (mul-terms L1 L2)
    (if (empty-termlist? L1)
        (the-empty-termlist)
        (add-terms (mul-term-by-all-terms (term-list-order L1) (first-term L1) L2)
                   (mul-terms (rest-terms L1) L2))))
  (define (mul-term-by-all-terms o1 c1 L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let* ((t2 (first-term L))
               (new-order (+ o1 (term-list-order L))))
          (adjoin-term
           new-order
           (mul c1 t2)
           (mul-term-by-all-terms o1 c1 (rest-terms L))))))

  ;; zero test
  (define (=zero-all-terms? L)
    (cond ((empty-termlist? L) #t)
          ((not (=zero? (first-term L))) #f)
          (else (=zero-all-terms? (rest-terms L)))))
        
  ;; Negation
  (define (negate-terms L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let ((term (first-term L)))
          (adjoin-term (term-list-order L)
                       (negate term)
                       (negate-terms (rest-terms L))))))
  
  ;; Equality
  (define (equ-terms? L1 L2)
    (cond ((empty-termlist? L1) (empty-termlist? L2))
          ((empty-termlist? L2) #f)
          (else (and (equ? (first-term L1) (first-term L2))
                     (equ-terms? (rest-terms L1) (rest-terms L2))))))
  
  ;; Constant
  (define (get-dense-constant L)
    (cond ((empty-termlist? L) zero)
          ((= (term-list-order L) 0) (first-term L))
          (else (get-dense-constant (rest-terms L)))))
  
  ;; Coercion
  (define (dense-terms->sparse-terms L)
    ((get 'make-from-coeffs 'sparse-terms) L))
  
  (define (count-zero-terms L)
    (if (empty-termlist? L)
        0
        (+ (if (=zero? (first-term L)) 1 0)
           (count-zero-terms (rest-terms L)))))
  
  (define (to-best-representation L)
    (if (store-as-sparse? (term-list-order L) (count-zero-terms L))
        (dense-terms->sparse-terms L)
        (tag L)))
  
  ;; interface to rest of the system
  (define (tag t) (attach-tag 'dense-terms (ensure-valid-term-list t)))
  (put 'add '(dense-terms dense-terms) 
       (lambda (t1 t2) (to-best-representation (add-terms t1 t2))))
  (put 'mul '(dense-terms dense-terms) 
       (lambda (t1 t2) (to-best-representation (mul-terms t1 t2))))
  (put 'equ? '(dense-terms dense-terms) equ-terms?)
  (put '=zero? '(dense-terms) =zero-all-terms?)
  (put 'negate '(dense-terms)
       (lambda (t) (to-best-representation (negate-terms t))))
  (put 'get-constant '(dense-terms) get-dense-constant)
  (put 'make-from-terms 'dense-terms
       (lambda (terms) (tag (make-from-terms terms))))
  (put 'make-from-coeffs 'dense-terms
       (lambda (coeffs) (tag (make-from-coeffs coeffs))))
  (put-coercion 'dense-terms 'sparse-terms dense-terms->sparse-terms)
  'done)

(define (get-constant term-list)
  (apply-generic 'get-constant term-list))
Finally we have the polynomial package that wraps the sparse and dense term-list representations. This basically takes the polynomial operations from the preceding exercises, strips out the term-list operations and replaces them with calls to the appropriate generic operations (i.e. add-terms becomes add, mul-terms becomes mul, and so on. And then we add a few tweaks...
  • There's the previously noted changes to same-variable? and the addition of select-variable to support zero-order "unbound" polynomials
  • We add an implementation of equ? to support coercion from polynomial to complex.
  • We add the coercion operation itself. Note that this takes the coefficient of the order 0 term from the underlying term-list and raises it to a complex so that it's of the correct type. Note that the coefficient may already be a complex (or may in fact be a polynomial itself), so we only try to raise-to a complex if it's lower than a complex.
One thing I didn't do in the polynomial package was to automatically select the term-list representation on construction. I decided to allow the clients of the system to determine the initial representation for their polynomials (or I was being lazy - you decide!).
;;;
;;; Polynomial wrapper package
;;; 
(define (install-polynomial-package)
  ;; internal procedures
  ;; representation of poly
  (define (make-poly variable term-list)
    (cons variable term-list))
  (define (make-from-coeffs variable coeffs)
    (make-poly variable
               ((get 'make-from-coeffs 'dense-terms) coeffs)))
  (define (make-from-terms variable terms)
    (make-poly variable
               ((get 'make-from-terms 'sparse-terms) terms)))
  
  (define (variable p) (car p))
  (define (term-list p) (cdr p))
  
  ;; variable tests and selection
  (define (variable? x) (symbol? x))
  (define (same-variable? v1 v2)
    (and (variable? v1)
         (variable? v2)
         (or (eq? v1 v2)
             (eq? v1 'unbound)
             (eq? v2 'unbound))))

  (define (select-variable p1 p2)
    (let ((v1 (variable p1)))
      (if (eq? v1 'unbound)
          (variable p2)
          v1)))

  ;; procedures used by add-poly
  (define (add-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (select-variable p1 p2)
                   (add (term-list p1)
                        (term-list p2)))
        (error "Polys not in same var -- ADD-POLY"
               (list p1 p2))))

  ;; procedures used by mul-poly
  (define (mul-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (select-variable p1 p2)
                   (mul (term-list p1)
                        (term-list p2)))
        (error "Polys not in same var -- MUL-POLY"
               (list p1 p2))))

  ;; Subtraction
  (define (sub-poly p1 p2)
    (add-poly p1 (negate-poly p2)))
  
  ;; zero test
  (define (=zero-poly? p)
    (=zero? (term-list p)))
        
  ;; Negation
  (define (negate-poly p)
    (make-poly (variable p)
               (negate (term-list p))))
  
  ;; Equality
  (define (equ-poly? p1 p2)
    (and (same-variable? (variable p1) (variable p2))
         (equ? (term-list p1) (term-list p2))))
  
  ;; Coercion
  (define (polynomial->complex p)
    (let ((constant (get-constant (term-list p))))
      (if (is-lower? constant 'complex)
          (raise-to 'complex constant)
          constant)))

  ;; interface to rest of the system
  (define (tag p) (attach-tag 'polynomial p))
  (put 'add '(polynomial polynomial) 
       (lambda (p1 p2) (tag (add-poly p1 p2))))
  (put 'mul '(polynomial polynomial) 
       (lambda (p1 p2) (tag (mul-poly p1 p2))))
  (put 'equ? '(polynomial polynomial) equ-poly?)
  (put '=zero? '(polynomial) =zero-poly?)
  (put 'negate '(polynomial)
       (lambda (p) (tag (negate-poly p))))
  (put 'sub '(polynomial polynomial)
       (lambda (p1 p2) (tag (sub-poly p1 p2))))
  (put 'make 'polynomial
       (lambda (var terms) (tag (make-poly var terms))))
  (put 'make-from-terms 'polynomial
       (lambda (variable terms) (tag (make-from-terms variable terms))))
  (put 'make-from-coeffs 'polynomial
       (lambda (variable coeffs) (tag (make-from-coeffs variable coeffs))))
  (put-coercion 'polynomial 'complex polynomial->complex)
  'done)

(define (make-polynomial-from-coeffs variable coeffs)
  ((get 'make-from-coeffs 'polynomial) variable coeffs))

(define (make-polynomial-from-terms variable terms)
  ((get 'make-from-terms 'polynomial) variable terms))
We have a special constructor for zero-order polynomials. As noted in the previous discussion, we're using a special marker, 'unbound, for zero-order polynomials that can participate in arithmetic operations with any other polynomial, regardless of their variable. This encapsulates the marker to support coercion from complex to polynomial types without complex having to know about 'unbound.
(define (make-zero-order-polynomial-from-coeff coeff)
  ((get 'make-from-coeffs 'polynomial) 'unbound (list coeff)))


;;;
;;; Install the Polynomial Packages
;;;
(install-term-package)
(install-sparse-terms-package)
(install-dense-terms-package)
(install-polynomial-package)
Okay, that's a big enough post for now... We'll take it for a spin in the next one!