2012-03-26

SICP Exercise 2.88: Subtracting Polynomials

Extend the polynomial system to include subtraction of polynomials. (Hint: You may find it helpful to define a generic negation operation.)
The authors suggest that defining a generic negation operation may be useful in solving this exercise. Why would this be?

Well, if we consider what's required in order to subtract one polynomial from another then all will become clear. Assuming we have polynomials p1 and p2 and we want to calculate p1 - p2 then, similar to the addition of polynomials, we want to perform this term-wise:
  • If there are terms of the same order in both p1 and p2 then we'll need a corresponding term in the resulting polynomial. The coefficient of this term will be the result of subtracting the coefficient of p2's term from the coefficient of p1's term.
  • If only p1 has a term of a particular order then we need an identical term in the resulting polynomial, so we can just add it into the resulting polynomial directly.
  • If only p2 has a term of a particular order then we need the negation of that term in the resulting polynomial.
So we could define an operation, sub-terms, that implements the steps above in a manner analogous to add-terms, and then define a procedure sub-poly which uses this in the same way as add-poly uses add-terms, using our generic negation operation to negate terms in p2 for which there is no corresponding term of the same order in p1.

Or we could be a bit smarter...

For any two values, a and b, that are of types representable in our system we know that a - b = a + (-b). In other words, subtracting b from a will give the same result as adding the negation of b to a.

When it comes to dealing with the terms of our polynomials this means that we can calculate the result of subtracting one term, t2, from another term, t1, by adding the result of negating t2 to t1. We can extend this to calculating the result of subtracting the polynomial p2 from p1 by negating all of the terms of p2, and then adding this new (negated) polynomial to p1. We already have an operation, add-poly, which will perform the addition. As a result, if we define a generic negation operation then implementing sub-poly simply becomes a case of evaluating add-poly with the first polynomial and the negation of the second polynomial.

So how do we define negation? Well, the generic operation itself is defined in the usual manner. We then just need to consider what to do for each of the different types:
  • integer: Given an integer i we can calculate its negation by evaluating (- i), and then all we need to do is tag the results of this evaluation.
  • rational: To negate a rational number we simply negate its numerator, which we can do in the same way as for integer, and then create a new rational number with the negated numerator and original denominator.
  • real: A real number can be negated in the same way as for an integer.
  • complex: To negate a complex number consider that, in order to calculate the negation of the complex number a + bi, we need to calculate -(a + bi). This is equivalent to (-a) + (-b)i. As a result we can calculate the negation by negating both the real and imaginary parts and constructing a new complex number from these negated components.
  • polynomial: Finally, to negate a polynomial, we simply iterate through its terms and negate each one in turn using the generic negation operation. Once we've done this we create a new polynomial with the same indeterminate as the original and with the negated term list.
Okay, given that, here's the changes required to provide a generic negation operation:
(define (install-integer-package)
  …
  (put 'negate '(integer)
       (lambda (x) (tag (- x))))
  …
  'done)

(define (install-rational-package)
  ;; internal procedures
  …
  (define (negate-rat x)
    (make-rat (- (numer x)) (denom x)))
  …
  ;; interface to rest of the system
  …
  (put 'negate '(rational)
       (lambda (x) (tag (negate-rat x))))
  …
  'done)

(define (install-real-package)
  …
  (put 'negate '(real)
       (lambda (x) (tag (- x))))
  …
  'done)

(define (install-complex-package)
  …
  ;; internal procedures
  …
  (define (negate-complex z)
    (make-from-real-imag (negate (complex-real-part z))
                         (negate (complex-imag-part z))))
  …
  ;; interface to rest of the system
  …
  (put 'negate '(complex)
       (lambda (z) (tag (negate-complex z))))
  …
  'done)

(define (install-polynomial-package)
  ;; internal procedures
  …
  ;; 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))))))
  (define (negate-poly p)
    (make-poly (variable p)
               (negate-terms (term-list p))))
  …
  ;; interface to rest of the system
  …
  (put 'negate '(polynomial)
       (lambda (p) (tag (negate-poly p))))
  …
  'done)
Here's the negation in action:
> (negate (make-integer 42))
(integer . -42)
> (negate (make-rational 13 27))
(rational -13 . 27)
> (negate (make-real 3.14159))
(rational -3537115888337719 . 1125899906842624)
> (negate (make-complex-from-real-imag (make-real 3.5) (make-integer -32)))
(complex rectangular (rational -7 . 2) integer . 32)
> (define inner-poly (make-polynomial 
                      'y
                      (list (make-polynomial-term 2 (make-integer 1))
                            (make-polynomial-term 1 (make-rational -3 2))
                            (make-polynomial-term 0 (make-integer 42)))))
> (define outer-poly (make-polynomial
                      'x
                      (list (make-polynomial-term 4 (make-rational 7 2))
                            (make-polynomial-term 3 (make-complex-from-real-imag
                                                     (make-integer 42)
                                                     (make-integer -13)))
                            (make-polynomial-term 2 inner-poly)
                            (make-polynomial-term 0 (make-integer 5)))))
> (negate outer-poly)
(polynomial x (4 (rational -7 . 2))
              (3 (complex rectangular (integer . -42) integer . 13))
              (2 (polynomial y (2 (integer . -1))
                               (1 (rational 3 . 2))
                               (0 (integer . -42))))
              (0 (integer . -5)))
With that working, we just implement sub-poly as noted above:
(define (install-polynomial-package)
  …
  ;; Subtraction
  (define (sub-poly p1 p2)
    (add-poly p1 (negate-poly p2)))
  …
  ;; interface to rest of the system
  …
  (put 'sub '(polynomial polynomial)
       (lambda (p1 p2) (tag (sub-poly p1 p2))))
  …
  'done)
Now we can subtract one polynomial from another... Subtracting 5x4 + 4x3 + x + 1 from 4x4 + 2x2 + 1 should give -x4 - 4x3 + 2x2 - x:
> (sub (make-polynomial 'x
                        (list (make-polynomial-term 4 (make-integer 4))
                              (make-polynomial-term 2 (make-integer 2))
                              (make-polynomial-term 0 (make-integer 1))))
       (make-polynomial 'x
                        (list (make-polynomial-term 4 (make-integer 5))
                              (make-polynomial-term 3 (make-integer 4))
                              (make-polynomial-term 1 (make-integer 1))
                              (make-polynomial-term 0 (make-integer 1)))))
(polynomial x (4 (integer . -1))
              (3 (integer . -4))
              (2 (integer . 2))
              (1 (integer . -1)))
Looks good!

2012-03-24

SICP Exercise 2.87 Addendum: Making Polynomials

If you were paying close attention at the end of exercise 2.87 you would have noticed that I created polynomials for testing =zero? as follows:
(make-polynomial 'x '())
(make-polynomial 'x (list (list 4 (make-integer 3))
                          (list 2 (make-integer 1))
                          (list 0 (make-real 2.3))))
(make-polynomial 'x (list (list 3 (make-real 0))
                          (list 2 (make-rational 0 4))
                          (list 1 (make-integer 0))))
Do you notice anything about the term list?

Yep, I had to create the term list representation directly as the operations within the polynomial package for doing so are not exposed externally. This means that I, and anyone using the system, needs to have knowledge of the internal representation used for terms and term lists in order to create polynomials. They would need to know that a term is a list consisting of the order of the indeterminate that the term is for and the coefficient to apply. They would also need to know that a term list is a list of such terms, ordered in decreasing value of the order of the terms.

Exposing the internal representation in this manner has a further implication. It means that once our system is in use we will be unable to change the internal representation without having to find all clients of the system and update them to use the new representation.

Yuck.

We can address these issues by encapsulating the internal representation of terms. To do this we'll need a way for clients to create valid term representations for passing to make-polynomial. We'll also want to update make-polynomial so that clients can pass us in any list of terms and we'll reorganize them into our desired internal representation.

The first of these steps is straightforward. The polynomial package already has an internal operation, make-term for constructing polynomial terms. We just need to expose this externally in the normal manner.

As for accepting any list of polynomials and converting them into a valid representation, let's consider what we'll have to do. We'll need to go through the passed-in list of terms and build up the valid internal representation as we go. To do this we'll need to start with a valid internal representation of an empty term list (which we can get by evaluating (the-empty-termlist)). We then need to iterate through the passed-in terms and, for each one, insert it into the correct position in the term list we're building up. We should also cope with the situation where there are two or more terms passed in with the same order. We could either raise an error in this case, or add the terms together. I chose the latter.

Here's the updates required to do all this:
(define (install-polynomial-package)
  …
  (define (make-poly variable term-list)
    (cons variable (build-term-list term-list (the-empty-termlist))))
  …
  (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-term-list terms result)
    (if (null? terms)
        result
        (build-term-list (cdr terms) (insert-term (car terms) result))))
  …
  (put 'make 'polynomial-term make-term)
  …
  'done)

(define (make-polynomial-term order coeff)
  ((get 'make 'polynomial-term) order coeff))
We can then make polynomials as follows:
> (make-polynomial 'x
                   (list (make-polynomial-term 4 (make-integer 3))
                         (make-polynomial-term 2 (make-integer 1))
                         (make-polynomial-term 0 (make-real 2.3))))
(polynomial x (4 (integer . 3)) (2 (integer . 1)) (0 (real . 2.3)))
> (make-polynomial 'y
                   (list (make-polynomial-term 2 (make-integer 1))
                         (make-polynomial-term 0 (make-real 2.3))
                         (make-polynomial-term 4 (make-integer 3))))
(polynomial y (4 (integer . 3)) (2 (integer . 1)) (0 (real . 2.3)))
> (make-polynomial 'z
                   (list (make-polynomial-term 0 (make-integer 0))
                         (make-polynomial-term 3 (make-real 0.0))
                         (make-polynomial-term 2 (make-rational 1 2))))
(polynomial z (2 (rational 1 . 2)))
> (make-polynomial 'a
                   (list (make-polynomial-term 2 (make-integer 1))
                         (make-polynomial-term 3 (make-real 3.5))
                         (make-polynomial-term 2 (make-rational 1 2))))
(polynomial a (3 (real . 3.5)) (2 (rational 3 . 2)))

2012-03-23

SICP Exercise 2.87: Testing Zero for Polynomials

Install =zero? for polynomials in the generic arithmetic package. This will allow adjoin-term to work for polynomials with coefficients that are themselves polynomials.

The polynomial Package

In section 2.5.3 the authors provide us with various operations for manipulating polynomials and their terms, along with the skeleton of a polynomial package. Before we can start into the exercise we need to put all this together so that we have a package we can install into the system we've been building from section 2.5.1.

To save you the hassle of doing it yourself, here's the basic polynomial package:
(define (install-polynomial-package)
  ;; internal procedures
  ;; representation of poly
  (define (make-poly variable term-list)
    (cons variable 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 (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))
  (define (make-term order coeff) (list order coeff))
  (define (order term) (car term))
  (define (coeff term) (cadr term))

  ;; 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)))
             (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-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 (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))))))
  
  ;; 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 'make 'polynomial
       (lambda (var terms) (tag (make-poly var terms))))
  'done)

(define (make-polynomial var terms)
  ((get 'make 'polynomial) var terms))
Don't forget to install this by evaluating (install-polynomial-package) of course, or you won't get very far!

Zero Testing for Polynomials

Now onto the exercise itself!

In order to implement =zero? for the polynomial package we need to ask ourselves how we can test whether a polynomial is zero or not. We can't determine equality with zero by looking at the indeterminate itself. However, we can determine it from the coefficients.

If there are no coefficients then the polynomial is equal to zero. On the other hand, when there are coefficients, we can only say that the polynomial is equal to zero iff all of the coefficients are equal to zero. As soon as one of them is non-zero we're unable to tell whether or not the polynomial is equal to zero until we bind the indeterminate to a particular value. Of course, as the coefficients are types from our generic arithmetic system we can test these for equality with zero by recursively evaluating =zero? on each of them.

So we have a recursive approach. If the list of coefficients is empty then the polynomial is equal to zero. Otherwise, we test whether or not the first coefficient is equal to zero. If it is then the polynomial may be equal to zero, but we need to test the remaining coefficients for equality with zero before we can confirm this; if it isn't then the polynomial isn't equal to zero.

Let's put this into code:
(define (install-polynomial-package)
  …
  (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)))
  …
  (put '=zero? '(polynomial) =zero-poly?)
  …
  'done)
...and let's give it a spin:
> (=zero? (make-polynomial 'x '()))
#t
> (=zero? (make-polynomial 'x (list (list 4 (make-integer 3))
                                    (list 2 (make-integer 1))
                                    (list 0 (make-real 2.3)))))
#f
> (=zero? (make-polynomial 'x (list (list 3 (make-real 0))
                                    (list 2 (make-rational 0 4))
                                    (list 1 (make-integer 0)))))
#t