2012-04-30

SICP Exercise 2.89 Redux: Stripping Out the Deadwood

Define procedures that implement the term-list representation described above as appropriate for dense polynomials.
Right at the start of exercise 2.89 I suggested that, in anticipation of exercise 2.90 it might be best to produce the dense polynomial package by making as little change to the general structure of the polynomial package as possible and that this would mean leaving in some procedures that could be stripped out.

Well, I've just started into the exercise (yes, I know, it's been over a month since I posted the last one, but I've been busy!) and come to the conclusion that this won't really help much due to the extent of the changes that'll be required. So here's exercise 2.89 with the "procedures that could be stripped out" actually stripped out:
(define (install-polynomial-package)
  ;; internal procedures
  ;; representation of poly
  (define (make-poly variable term-list)
    (cons variable (strip-leading-zeros term-list)))
  (define (variable p) (car p))
  (define (term-list p) (cdr p))
  
  ;; procedures same-variable? and variable? from section 2.3.2
  (define (variable? x) (symbol? x))
  (define (same-variable? v1 v2)
    (and (variable? v1) (variable? v2) (eq? v1 v2)))

  ;; representation of terms and term lists
  (define zero-term (make-integer 0))
  (define (adjoin-term term-order term-coeff term-list)
    (cond ((=zero? term-coeff) term-list)
          ((= term-order (+ 1 (order term-list)))
           (cons term-coeff term-list))
          ((> term-order (order term-list))
           (adjoin-term term-order
                        term-coeff
                        (cons zero-term 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 (order term-list)
    (- (length term-list) 1))

  (define (strip-leading-zeros term-list)
    (cond ((empty-termlist? term-list) (the-empty-termlist))
          ((not (=zero? (first-term term-list))) term-list)
          (else (strip-leading-zeros (rest-terms term-list)))))
  
  ;; procedures used by add-poly
  (define (add-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (variable p1)
                   (add-terms (term-list p1)
                              (term-list p2)))
        (error "Polys not in same var -- ADD-POLY"
               (list p1 p2))))
  (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 (order L1))
                  (o2 (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-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (variable p1)
                   (mul-terms (term-list p1)
                              (term-list p2)))
        (error "Polys not in same var -- MUL-POLY"
               (list p1 p2))))
  (define (mul-terms L1 L2)
    (if (empty-termlist? L1)
        (the-empty-termlist)
        (add-terms (mul-term-by-all-terms (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 (order L))))
          (adjoin-term
           new-order
           (mul c1 t2)
           (mul-term-by-all-terms o1 c1 (rest-terms L))))))

  ;; Subtraction
  (define (sub-poly p1 p2)
    (add-poly p1 (negate-poly p2)))
  
  ;; 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)))))
  (define (=zero-poly? p)
    (=zero-all-terms? (term-list p)))
        
  ;; Negation
  (define (negate-terms L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let ((term (first-term L)))
          (adjoin-term (order L)
                       (negate term)
                       (negate-terms (rest-terms L))))))
  (define (negate-poly p)
    (make-poly (variable p)
               (negate-terms (term-list p))))
  
  ;; 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 '=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))))
  'done)
Right, now back to exercise 2.90...

2012-04-08

SICP Exercise 2.89: Representing Dense Polynomials

Define procedures that implement the term-list representation described above as appropriate for dense polynomials.
In anticipation of the next exercise, let's try to do this with as little change to the general structure of the polynomial package as possible. That way it should be less of a step to supporting two different representations of term lists when we get there. This means we're going to leave in some procedures that could be stripped out (I'm looking at you, make-term and coeff). We will, of course, have to make a few changes to the existing system, but we'll try to do these in such a way that we can support both dense and sparse polynomials easily after the changes.

Dense Polynomial Term Manipulation

First let's deal with the term operations. In our sparse polynomial representation each term consisted of both the order of the term and the coefficient. However, in the discussion on the representation of term lists for dense polynomials presented in the book the authors note that dense polynomial term lists "are most efficiently represented as lists of the coefficients", with coefficients of zero used to indicate zero terms.

As the terms for dense polynomials are just represented as the coefficients we could remove the operations make-poly and coeff. However, remembering that we're trying to minimise the set of changes we produce, let's retain these procedures with the same operands. For dense polynomials make-poly will just ignore the order operand and return the coeff operand unchanged, while coeff will just return the term operand unchanged:
(define (make-term order coeff) coeff)
(define (coeff term) term)
The remaining operation for manipulating terms, order, runs into a spot of bother with the dense polynomial representation. As the representation of a term is just its coefficient it's not possible to determine the order of a term from the term itself. However, the authors note that "the order of a term in this representation is the length of the sublist beginning with that term's coefficient, decremented by 1." So we can still implement an order procedure for dense polynomials, but it will need to work on the full term list, not just a single term. It can be implemented as follows:
(define (order termlist)
  (- (length term-list) 1))

Dense Polynomial Term List Manipulation

Now let's move on to manipulating term lists for dense polynomials. We've already noted that the authors suggest representing the term lists as list of the coefficients. So, just like sparse polynomial term lists, dense polynomial term lists can be represented as Scheme lists. Of course the actual representation of the terms differ, but the data structure used to hold the terms is unchanged.

Note that, in order that we can determine which coefficient applies to which order of the indeterminate this representation requires that the coefficients be ordered in a predictable manner, and the text implies that we use an identical ordering as is used for the sparse polynomial term lists. I.e. from highest order term to lowest order term.

Given all this, we can create empty term lists, test for empty term lists and obtain the first and rest of the terms in the same manner as we do for sparse polynomials and so the-empty-termlist, empty-termlist?, first-term and rest-terms remain unchanged.

When it comes to adding terms to term lists, however, a bit more work is required. We need to ensure that the new term being adjoined to a term list is one order higher than the current order of the term list, and zero-pad the head of the list if it isn't. As we can't determine the order of the term from the term alone we'll now have to pass the order to adjoin-term. We'll then have to check the order of the term list and, if it's not one less than the new term's order, add a zero term to the head of the list and try again. And, because I'm paranoid, I'm also going to add a check that the term we're trying to adjoin is definitely of a higher order than the current term list...

Note that, instead of passing the order of the term and the term itself, we'll pass the order and coefficient of the term as the two operands. We can then create the term within adjoin-term using make-term, rather than creating it externally and passing it in. That makes little difference to dense polynomials, as the term is the coefficient. For sparse polynomials it makes the interface feel a little cleaner, as we won't be passing in both the order and the term, which contains the order.

Anyway, here's adjoin-term for dense polynomials:
(define zero-term (make-integer 0))
(define (adjoin-term term-order term-coeff term-list)
  (cond ((=zero? term-coeff) term-list)
        ((= term-order (+ 1 (order term-list)))
         (cons (make-term term-order term-coeff) term-list))
        ((> term-order (order term-list))
         (adjoin-term term-order
                      term-coeff
                      (cons zero-term term-list)))
        (else (error "Cannot adjoin term of lower order than term list -- ADJOIN-TERM"
                     (list term-order term-coeff term-list)))))

Updating the Dense Polynomial Operations

We've done our best above to minimise the changes we make to the polynomial package. However, we still had to modify order so that it takes the term list instead of just a single term and we've split the term operand for adjoin-term into both the order and the coefficient. This of course means that we need to modify the usages of these operations appropriately.

For adding dense polynomials we can update add-terms so that it calculates the orders of the terms at the same time it's pulling out the terms at the heads of each list. We'll use let* to do this succinctly. Then we can update the calls to adjoin-term to pass through the order:
(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 (order L1))
                (o2 (order L2)))
           (cond ((> o1 o2)
                  (adjoin-term
                   o1 (coeff t1) (add-terms (rest-terms L1) L2)))
                 ((< o1 o2)
                  (adjoin-term
                   o2 (coeff t2) (add-terms L1 (rest-terms L2))))
                 (else
                  (adjoin-term
                   o1
                   (add (coeff t1) (coeff t2))
                   (add-terms (rest-terms L1)
                              (rest-terms L2)))))))))
Multiplying polynomials requires a more complex process. We go through each term in the first polynomial's term list and, for each one, produce a new term list that's the result of multiplying each term in the second polynomial's term list by this term. We add up all these term lists to produce the term list for the polynomial that results from the multiplication. The outer loop for this is implemented by mul-terms, while the inner loop (i.e. the multiplication of each term in the second polynomial's term list by a term in the first polynomial is implemented by mul-term-by-all-terms.

Looking at the inner loop first, we can see that this needs to know the order of the term from the first polynomial in order to calculate the results of multiplying this term by the first term in the second polynomial's term list. Like adjoin-term this means that we will need to pass in the order of the term as it's no longer possible to calculate it from the term itself. In fact we can apply a change to mul-term-by-all-terms that's similar to the one we applied to adjoin-term. I.e. pass in the order and the coefficient of the term as separate operands.

This then has a knock-on effect on mul-terms, where we need to pass these in to the call to mul-term-by-all-terms.

Here's the updated operations:
(define (mul-terms L1 L2)
  (if (empty-termlist? L1)
      (the-empty-termlist)
      (add-terms (mul-term-by-all-terms (order L1) (coeff (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 (order L))))
        (adjoin-term
         new-order
         (mul c1 (coeff t2))
         (mul-term-by-all-terms o1 c1 (rest-terms L))))))
There's one more operation to consider, sub-poly. However, this just calls through to add-poly with the first polynomial and the negation of the second polynomial that it was invoked with. So sub-poly remains unchanged. We need to make a minor change to negate-terms to account for the change in the order and adjoin-term operations though:
(define (negate-terms L)
  (if (empty-termlist? L)
      (the-empty-termlist)
      (let ((term (first-term L)))
        (adjoin-term (order L)
                     (negate (coeff term))
                     (negate-terms (rest-terms L))))))

Creating Dense Polynomials

The last thing we need to deal with for dense polynomials is constructing them.

We should note that the list-of-coefficients representation allows us to use a simpler interface for constructing dense polynomials than for sparse polynomials. Provided we stipulate that the list of coefficients must be ordered from highest order coefficient to lowest, has a coefficient for order 0 and has zero coefficients for zero terms then all we require are the indeterminate and the list of the coefficients; we don't need the orders of the terms any more as this is implicit. As a result, we won't need to expose an operation for making terms for dense polynomials such as was added in my addendum to exercise 2.87.

In that addendum I also modified make-poly from the implementation provided in the book so that it reordered the terms passed to it to ensure that they were ordered from high to low order. Given that we're making the ordering of coefficients into part of the creation operation's contract we no longer need to perform this reordering. However, we'll still want to pre-process the provided set of coefficients to remove any leading zero coefficients. Leading zero coefficients do not contribute to the evaluation of the polynomial, and so can be safely ignored.

We can do this quite simply by iterating through the term list, checking the head of the list at each iteration with =zero?, and returning the list as soon as we encounter a non-zero entry. Obviously we'll also need to check to see if we've got an empty term list - as could happen if we're passed either an entirely empty term list or one in which all coefficients are zero.

Here's make-poly and the accompanying helper operation:
(define (make-poly variable term-list)
  (cons variable (strip-leading-zeros term-list)))
(define (strip-leading-zeros term-list)
  (cond ((empty-termlist? term-list) (the-empty-termlist))
        ((not (=zero? (first-term term-list))) term-list)
        (else (strip-leading-zeros (rest-terms term-list)))))
Note that, although the term list representations differ, the make-poly implementations for both dense and sparse polynomials create the same "outer" structure for representing polynomials: a pair with the indeterminate as the first item and the term list as the second item. This means that the variable and term-list operations can remain unchanged. Also, as we're not changing the representation of the indeterminate, variable? and same-variable? can also remain as is.

Putting It All Together

We've now gone through the entirety of the polynomial package and produced a minimal set of changes to convert it from sparse polynomial to dense polynomial representation. You could pull the dense polynomial version together from above, but to save you the hassle, here it its:
(define (install-polynomial-package)
  ;; internal procedures
  ;; representation of poly
  (define (make-poly variable term-list)
    (cons variable (strip-leading-zeros term-list)))
  (define (variable p) (car p))
  (define (term-list p) (cdr p))
  
  ;; procedures same-variable? and variable? from section 2.3.2
  (define (variable? x) (symbol? x))
  (define (same-variable? v1 v2)
    (and (variable? v1) (variable? v2) (eq? v1 v2)))

  ;; representation of terms and term lists
  (define zero-term (make-integer 0))
  (define (adjoin-term term-order term-coeff term-list)
    (cond ((= term-order (+ 1 (order term-list)))
           (cons (make-term term-order term-coeff) term-list))
          ((> term-order (order term-list))
           (adjoin-term term-order
                        term-coeff
                        (cons zero-term 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 (make-term order coeff) coeff)
  (define (order term-list)
    (- (length term-list) 1))
  (define (coeff term) term)

  (define (strip-leading-zeros term-list)
    (cond ((empty-termlist? term-list) (the-empty-termlist))
          ((not (=zero? (first-term term-list))) term-list)
          (else (strip-leading-zeros (rest-terms term-list)))))
  
  ;; procedures used by add-poly
  (define (add-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (variable p1)
                   (add-terms (term-list p1)
                              (term-list p2)))
        (error "Polys not in same var -- ADD-POLY"
               (list p1 p2))))
  (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 (order L1))
                  (o2 (order L2)))
             (cond ((> o1 o2)
                    (adjoin-term
                     o1 (coeff t1) (add-terms (rest-terms L1) L2)))
                   ((< o1 o2)
                    (adjoin-term
                     o2 (coeff t2) (add-terms L1 (rest-terms L2))))
                   (else
                    (adjoin-term
                     o1
                     (add (coeff t1) (coeff t2))
                     (add-terms (rest-terms L1)
                                (rest-terms L2)))))))))

  ;; procedures used by mul-poly
  (define (mul-poly p1 p2)
    (if (same-variable? (variable p1) (variable p2))
        (make-poly (variable p1)
                   (mul-terms (term-list p1)
                              (term-list p2)))
        (error "Polys not in same var -- MUL-POLY"
               (list p1 p2))))
  (define (mul-terms L1 L2)
    (if (empty-termlist? L1)
        (the-empty-termlist)
        (add-terms (mul-term-by-all-terms (order L1) (coeff (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 (order L))))
          (adjoin-term
           new-order
           (mul c1 (coeff t2))
           (mul-term-by-all-terms o1 c1 (rest-terms L))))))

  ;; Subtraction
  (define (sub-poly p1 p2)
    (add-poly p1 (negate-poly p2)))
  
  ;; 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)))))
  (define (=zero-poly? p)
    (=zero-all-terms? (term-list p)))
        
  ;; Negation
  (define (negate-terms L)
    (if (empty-termlist? L)
        (the-empty-termlist)
        (let ((term (first-term L)))
          (adjoin-term (order L)
                       (negate (coeff term))
                       (negate-terms (rest-terms L))))))
  (define (negate-poly p)
    (make-poly (variable p)
               (negate-terms (term-list p))))
  
  ;; 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 '=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))))
  'done)

(define (make-polynomial var terms)
  ((get 'make 'polynomial) var terms))
Now we've got our dense polynomial package, let's give it a spin:
> (make-polynomial 'x '())
(polynomial x)
> (make-polynomial 'x (list (make-integer 3)
                            (make-integer 1)
                            (make-real 2.3)))
(polynomial x (integer . 3) (integer . 1) (real . 2.3))
> (make-polynomial 'x (list (make-real 0)
                            (make-rational 0 4)
                            (make-integer 0)
                            (make-complex-from-real-imag
                             (make-integer 0)
                             (make-integer 0))))
(polynomial x)
> (define inner-poly (make-polynomial 
                      'y
                      (list (make-integer 1)
                            (make-rational -3 2)
                            (make-integer 42))))
> (define outer-poly (make-polynomial
                      'x
                      (list (make-rational 7 2)
                            (make-complex-from-real-imag
                             (make-integer 42)
                             (make-integer -13))
                            inner-poly
                            (make-integer 0)
                            (make-integer 5))))
> (negate outer-poly)
(polynomial x (rational -7 . 2)
              (complex rectangular (integer . -42) integer . 13)
              (polynomial y (integer . -1) (rational 3 . 2) (integer . -42))
              (integer . 0)
              (integer . -5))
> (add (make-polynomial 'x
                        (list (make-integer 4)
                              (make-integer 0)
                              (make-integer 2)
                              (make-integer 0)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 5)
                              (make-integer 4)
                              (make-integer 0)
                              (make-integer 1)
                              (make-integer 1))))
(polynomial x (integer . 9)
              (integer . 4)
              (integer . 2)
              (integer . 1)
              (integer . 2))
> (add (make-polynomial 'x
                        (list (make-integer 5)
                              (make-integer 4)
                              (make-integer 3)
                              (make-integer 2)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 3)
                              (make-integer 2)
                              (make-integer 1))))
(polynomial x (integer . 5)
              (integer . 4)
              (integer . 6)
              (integer . 4)
              (integer . 2))
> (add (make-polynomial 'x
                        (list (make-integer 3)
                              (make-integer 2)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 5)
                              (make-integer 4)
                              (make-integer 3)
                              (make-integer 2)
                              (make-integer 1))))
(polynomial x (integer . 5)
              (integer . 4)
              (integer . 6)
              (integer . 4)
              (integer . 2))
> (sub (make-polynomial 'x
                        (list (make-integer 4)
                              (make-integer 0)
                              (make-integer 2)
                              (make-integer 0)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 5)
                              (make-integer 4)
                              (make-integer 0)
                              (make-integer 1)
                              (make-integer 1))))
(polynomial x (integer . -1)
              (integer . -4)
              (integer . 2)
              (integer . -1)
              (integer . 0))
> (sub (make-polynomial 'x
                        (list (make-integer 3)
                              (make-integer 2)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 5)
                              (make-integer 4)
                              (make-integer 3)
                              (make-integer 2)
                              (make-integer 1))))
(polynomial x (integer . -5)
              (integer . -4)
              (integer . 0)
              (integer . 0)
              (integer . 0))
> (mul (make-polynomial 'x
                        (list (make-integer 3)
                              (make-integer 2)
                              (make-integer 1)))
       (make-polynomial 'x
                        (list (make-integer 6)
                              (make-integer 5)
                              (make-integer 4)
                              (make-integer 3)
                              (make-integer 2)
                              (make-integer 1))))
(polynomial x (integer . 18)
              (integer . 27)
              (integer . 28)
              (integer . 22)
              (integer . 16)
              (integer . 10)
              (integer . 4)
              (integer . 1))