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 bynew-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 completediv-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 ofdiv-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
(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)))
No comments:
Post a Comment