Implement 2.93

This took way longer than expected. First, I noticed that mixed up
sparse and dense terms so I fixed that. Then I extended the rational
package to use generic operations. However, I did not work and it took
me two hours to find the root cause. The problem was that appply-generic
tried to raise raise and drop the rationals with polynomials. The
rational to scheme and rational to real procedures do not support
polynomials which caused the error. I changed these procedures so that
they only attempt coercion of the arguments are numbers. This is not
fully correct as it destroys our ability to have rationals of rationals,
for example, but good enough for our purposes.
main
Felix Martin 2020-12-05 12:48:56 -05:00
parent 47e979a30c
commit d25e14aaa2
1 changed files with 125 additions and 89 deletions

View File

@ -59,6 +59,10 @@
(define (tag x) x)
(define (scheme->rational x)
(make-rational x 1))
(define (gcd-scheme a b)
(if (and (integer? a) (integer? b))
(if (= b 0) (abs a) (gcd-scheme b (remainder a b)))
'nogcd))
(put 'add '(scheme-number scheme-number)
(lambda (x y) (tag (+ x y))))
(put 'sub '(scheme-number scheme-number)
@ -81,6 +85,7 @@
(lambda (x y) (atan x y)))
(put 'square-root '(scheme-number) sqrt)
(put 'raise 'scheme-number scheme->rational)
(put 'gcd '(scheme-number scheme-number) gcd-scheme)
(display "[install-scheme-number-package]\n")
'done)
@ -89,33 +94,44 @@
(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)))
(cons n d)))
(let ((g (gcd n d)))
(if (eq? g 'nogcd)
(cons n d)
(cons (/ n g) (/ d g)))))
(define (add-rat x y)
(make-rat (+ (* (numer x) (denom y))
(* (numer y) (denom x)))
(* (denom x) (denom y))))
(let ((new-n (add (mul (numer x) (denom y))
(mul (numer y) (denom x))))
(new-d (mul (denom x) (denom y))))
(make-rat new-n new-d)))
(define (sub-rat x y)
(make-rat (- (* (numer x) (denom y))
(* (numer y) (denom x)))
(* (denom x) (denom y))))
(make-rat (sub (mul (numer x) (denom y))
(mul (numer y) (denom x)))
(mul (denom x) (denom y))))
(define (mul-rat x y)
(make-rat (* (numer x) (numer y))
(* (denom x) (denom y))))
(make-rat (mul (numer x) (numer y))
(mul (denom x) (denom y))))
(define (div-rat x y)
(make-rat (* (numer x) (denom y))
(* (denom x) (numer y))))
(make-rat (mul (numer x) (denom y))
(mul (denom x) (numer y))))
(define (add3-rat x y z)
(add-rat (add-rat x y) z))
(define (equ? x y)
(= (* (numer x) (denom y))
(* (numer y) (denom x))))
(define (equ?-rat x y)
(equ? (mul (numer x) (denom y))
(mul (numer y) (denom x))))
(define (rational->real x)
(make-real (/ (numer x) (denom x))))
(let ((n (numer x))
(d (denom x)))
(cond
((and (number? n) (number? d))
(make-real (/ (numer x) (denom x))))
(else 'invalid))))
(define (rational->scheme x)
(make-scheme-number (inexact->exact (round (/ (numer x) (denom x))))))
(let ((n (numer x)) (d (denom x)))
(cond
((and (number? n) (number? d))
(make-scheme-number (inexact->exact (round (/ (numer x) (denom x))))))
(else 'invalid))))
;; interface to rest of the system
(define (tag x) (attach-tag 'rational x))
(put 'add '(rational rational)
@ -128,7 +144,7 @@
(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 'equ? '(rational rational) equ?-rat)
(put '=zero? '(rational)
(lambda (x) (= (numer x) 0)))
(put 'negate '(rational)
@ -316,6 +332,10 @@
(define (arctan x y) (apply-generic 'arctan x y))
(define (square-root x) (apply-generic 'square-root x))
(define (negate x) (apply-generic 'negate x))
(define (gcd x y)
(if (procedure? (get 'gcd (list (type-tag x) (type-tag y))))
(apply-generic 'gcd x y)
'nogcd))
(install-scheme-number-package)
(install-rational-package)
@ -468,13 +488,15 @@
; All we have to do is update coerce-args to do consecutive raises
; to reach the target type.
(define (coerce-args target-type args)
; (display "COERCE-ARGS ") (display target-type) (display " ") (display args) (newline)
(define (coerce-arg arg)
; (display "COERCE-ARG ") (display arg) (newline)
(if (eq? (type-tag arg) target-type)
arg
(let ((raise (get 'raise (type-tag arg))))
(if (procedure? raise)
(raise (contents arg))
arg))))
(if (procedure? raise)
(raise (contents arg))
arg))))
(let ((coerced-args (map coerce-arg args)))
(if (equal? args coerced-args)
coerced-args ; no more raising possible
@ -501,14 +523,15 @@
; Implement drop to transform number to lowest possible representation
(define (drop x)
;(display "---------\ndrop ") (display x) (newline)
; (display "DROP ") (display x) (newline)
(if (has-tag? x)
(let ((project (get 'project (type-tag x))))
(if (procedure? project)
(let ((projected (project (contents x))))
(if (equ? projected x)
(drop projected)
x))
(cond
((eq? projected 'invalid) x)
((equ? projected x) (drop projected))
(else x)))
x))
x))
@ -517,7 +540,7 @@
;(assert (drop (make-complex-from-real-imag 3 0)) (make-scheme-number 3))
(define (apply-generic op . args)
;(display "-----\napply-generic ") (display op) (display " ") (display args) (newline)
; (display "APPLY-GENERIC ") (display op) (display " ") (display args) (newline)
(define (try-args args-list)
(if (null? args-list)
(error "No method for these types" (list op (map type-tag args)))
@ -527,6 +550,7 @@
(drop (apply proc args-contents))
(try-args (cdr args-list))))))
(define (coerce-to-arg arg)
; (display "COERCE-TO-ARG ") (display arg) (newline)
(coerce-args (type-tag arg) args))
(try-args (cons args (map coerce-to-arg args))))
@ -728,52 +752,50 @@
((get 'first-term (type-tag term-list)) (contents term-list)))
(define (adjoin-term term term-list)
; (display "ADJOIN-TERM ") (display term) (display term-list) (newline)
((get 'adjoin-term (type-tag term-list)) term (contents term-list)))
;; sparse implementations
(define (tag-sparse p) (attach-tag 'sparse p))
(define (make-poly-sparse variable term-list)
(cons variable (tag-sparse term-list)))
(define (first-term-sparse term-list)
;(display "first-term-sparse ") (display term-list) (newline)
(make-term (- (length term-list) 1)
(car term-list)))
(define (adjoin-term-sparse term term-list)
(let ((coeff-term (coeff term))
(order-term (order term))
(length-terms (length term-list)))
(cond
((= order-term length-terms) (tag-sparse (cons coeff-term term-list)))
((< order-term length-terms) (error "Cannot adjoin lower-order term to terms"))
(else (tag-sparse
(cons
coeff-term
(contents (adjoin-term-sparse
(make-term (- order-term 1) 0)
term-list))))))))
(put 'first-term 'sparse first-term-sparse)
(put 'adjoin-term 'sparse adjoin-term-sparse)
;; dense implementations
(define (tag-dense p) (attach-tag 'dense p))
(define (make-poly-dense variable term-list)
(cons variable (tag-dense term-list)))
(define (adjoin-term-dense term term-list)
(if (=zero? (coeff term))
(tag-dense term-list)
(tag-dense (cons term term-list))))
(define (first-term-dense term-list)
; (display "FIRST-TERM-DENSE ") (display term-list) (newline)
(make-term (- (length term-list) 1)
(car term-list)))
(define (first-term-dense term-list) (car term-list))
(define (adjoin-term-dense term term-list)
(let ((coeff-term (coeff term))
(order-term (order term))
(length-terms (length term-list)))
(cond
((= order-term length-terms) (tag-dense (cons coeff-term term-list)))
((< order-term length-terms) (error "Cannot adjoin lower-order term to terms"))
(else (tag-dense
(cons
coeff-term
(contents (adjoin-term-dense
(make-term (- order-term 1) 0)
term-list))))))))
(put 'first-term 'dense first-term-dense)
(put 'adjoin-term 'dense adjoin-term-dense)
(define (the-empty-termlist) (tag-sparse '()))
;; sparse implementations
(define (tag-sparse p) (attach-tag 'sparse p))
(define (make-poly-sparse variable term-list)
(cons variable (tag-sparse term-list)))
(define (adjoin-term-sparse term term-list)
(if (=zero? (coeff term))
(tag-sparse term-list)
(tag-sparse (cons term term-list))))
(define (first-term-sparse term-list) (car term-list))
(put 'first-term 'sparse first-term-sparse)
(put 'adjoin-term 'sparse adjoin-term-sparse)
(define (empty-termlist t)
(attach-tag (type-tag t) '()))
(define (rest-terms term-list)
@ -837,17 +859,17 @@
; (display "COERCE-TERM ") (display t) (newline)
(let ((c (coeff t))
(o (order t))
(new-poly (tag (make-poly-dense
(new-poly (tag (make-poly-sparse
source-var
(list (list (order t) 1))))))
; (display "NEW-POLY ") (display new-poly) (newline)
(if (eq? (type-tag c) 'polynomial)
(let ((new-poly (tag (make-poly-dense source-var (list (list o 1))))))
(let ((new-poly (tag (make-poly-sparse source-var (list (list o 1))))))
(scale-terms new-poly (contents (term-list c))))
(let ((sub-poly (tag (make-poly-dense
(let ((sub-poly (tag (make-poly-sparse
source-var
(list (list o c))))))
(cons 'dense (list (list 0 sub-poly)))))))
(cons 'sparse (list (list 0 sub-poly)))))))
(if (empty-termlist? terms)
terms
(add-terms (coerce-term (first-term terms))
@ -855,7 +877,7 @@
(define (coerce-poly p target-var)
; (display "COERCE-POLY ") (display p) (newline)
; (display "COERCE-POLY ") (display p) (display " TARGET ") (display target-var) (newline)
(if (eq? (variable p) target-var)
p
(let ((coercion-result (coerce-terms (term-list p) (variable p) target-var)))
@ -863,6 +885,7 @@
(make-poly target-var coercion-result))))
(define (coerce-polys p1 p2)
; (display "COERCE-POLYS ") (display p1) (display p2) (newline)
(let ((coercion-target-variable (get-coercion-target p1 p2)))
(if coercion-target-variable
(list (coerce-poly p1 coercion-target-variable)
@ -892,21 +915,25 @@
(list p1 p2))))
(define (mul-terms L1 L2)
; (display "MUL-TERMS") (display L1) (display L2) (newline)
(if (empty-termlist? L1)
L1
(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)
; (display "MUL-TERM-BY-ALL-TERMS ") (display t1) (display L) (newline)
(if (empty-termlist? L)
L
(let ((t2 (first-term L)))
; (display "T1 ") (display t1) (display " T2 ") (display t2) (newline)
(adjoin-term
(make-term (+ (order t1) (order t2))
(mul (coeff t1) (coeff t2)))
(mul-term-by-all-terms t1 (rest-terms L))))))
(define (mul-poly p1 p2)
; (display "MUL-POLY ") (display p1) (display p2) (newline)
(if (same-variable? (variable p1) (variable p2))
(make-poly (variable p1)
(mul-terms (term-list p1)
@ -986,52 +1013,61 @@
((get 'make 'poly-dense) var terms))
(define p1 (make-poly-sparse 'x (list 5 1)))
(define p1 (make-poly-dense 'x (list 5 1)))
(display p1) (newline)
(display (add p1 p1)) (newline)
(assert (add p1 p1) (make-poly-sparse 'x (list 10 2)))
(assert (add (make-poly-sparse 'x (list 2 2 0 1))
(make-poly-sparse 'x (list 1 2 3 2 3 6 6)))
(make-poly-sparse 'x (list 1 2 3 4 5 6 7)))
(assert (mul (make-poly-sparse 'x (list 1 1))
(make-poly-sparse 'x (list 1 1)))
(make-poly-sparse 'x (list 1 2 1)))
(assert (add p1 p1) (make-poly-dense 'x (list 10 2)))
(assert (add (make-poly-dense 'x (list 2 2 0 1))
(make-poly-dense 'x (list 1 2 3 2 3 6 6)))
(make-poly-dense 'x (list 1 2 3 4 5 6 7)))
(assert (mul (make-poly-dense 'x (list 1 1))
(make-poly-dense 'x (list 1 1)))
(make-poly-dense 'x (list 1 2 1)))
(display "\nex-2.90 - support sparse and dense polynomials\n")
(define p (make-poly-dense 'x '((100 2) (1 2))))
(define p (make-poly-sparse 'x '((100 2) (1 2))))
(display p) (newline)
(assert (add p p) (make-poly-dense 'x '((100 4) (1 4))))
(assert (add p p) (make-poly-sparse 'x '((100 4) (1 4))))
(assert (mul p p)
(make-poly-dense 'x '((200 4) (101 8) (2 4))))
(make-poly-sparse 'x '((200 4) (101 8) (2 4))))
(display "\nex-2.91 - divide\n")
(define p1 (make-poly-dense 'x '((5 1) (0 -1))))
(define p2 (make-poly-dense 'x '((2 1) (0 -1))))
(assert (mul p1 p1) (make-poly-dense 'x '((10 1) (5 -2) (0 1))))
(assert (mul p1 p2) (make-poly-dense 'x '((7 1) (5 -1) (2 -1) (0 1))))
(define p1 (make-poly-sparse 'x '((5 1) (0 -1))))
(define p2 (make-poly-sparse 'x '((2 1) (0 -1))))
(assert (mul p1 p1) (make-poly-sparse 'x '((10 1) (5 -2) (0 1))))
(assert (mul p1 p2) (make-poly-sparse 'x '((7 1) (5 -1) (2 -1) (0 1))))
(let ((result (div p1 p2)))
(assert (car result) (make-poly-dense 'x '((3 1) (1 1))))
(assert (cadr result) (make-poly-dense 'x '((1 1) (0 -1)))))
(assert (car result) (make-poly-sparse 'x '((3 1) (1 1))))
(assert (cadr result) (make-poly-sparse 'x '((1 1) (0 -1)))))
(display "\nex-2.92 - coerce polynomials\n")
(define p1 (make-poly-dense
(define p1 (make-poly-sparse
'x
(list (list 1 (make-poly-dense
(list (list 1 (make-poly-sparse
'y
'((2 1) (0 1)))))))
(define p2 (make-poly-dense
(define p2 (make-poly-sparse
'y
(list
(list 4 3)
(list 2 (make-poly-dense 'x '((2 1) (0 8)))))))
(list 2 (make-poly-sparse 'x '((2 1) (0 8)))))))
(display (add p1 p2)) (newline)
(display "\nex-2.93\n")
(display "\nex-2.93 - polynomial rationals\n")
(define p1 (make-poly-sparse 'x '((2 1)(0 1))))
(define p2 (make-poly-sparse 'x '((3 1)(0 1))))
(define rf (make-rational p2 p1))
(display rf) (newline)
(display (add rf rf)) (newline)
(display "\nex-2.94 - polynomial gcd\n")