home *** CD-ROM | disk | FTP | other *** search
/ ARM Club 3 / TheARMClub_PDCD3.iso / hensa / maths / b116_1 / jacal / elim < prev    next >
Text File  |  1993-10-24  |  9KB  |  285 lines
  1. ;;; JACAL: Symbolic Mathematics System.        -*-scheme-*-
  2. ;;; Copyright 1989, 1990, 1991, 1992, 1993 Aubrey Jaffer.
  3. ;;; See the file "COPYING" for terms applying to this program.
  4.  
  5. ;;;; Variable elimination
  6. (define (poly:elim poleqns vars)
  7.   (cond (math:trace
  8.      (display-diag 'eliminating:) (newline-diag)
  9.      (math:write (map var->expl vars) *output-grammar*)
  10.      (display-diag 'from:) (newline-diag)
  11.      (math:write (poleqns->licits poleqns) *output-grammar*)))
  12.   (do ((vs vars (cdr vs)) (polys poleqns) (poly #f))
  13.       ((null? vs)
  14.        (cond (math:trace (display-diag 'yielding:) (newline-diag)
  15.              (math:write polys *output-grammar*)))
  16.        polys)
  17.     (do ((var (car vs))
  18.      (pl polys (if (null? pl) 
  19.                (math:error 'not-enough-equations poleqns vars)
  20.                (cdr pl)))
  21.      (npl '() (cons (car pl) npl)))
  22.     ((poly:find-var? (car pl) var)
  23.      (set! poly (promote var (car pl)))
  24.      (do ((pls (cdr pl) (cdr pls)))
  25.          ((null? pls) (set! polys npl))
  26.        (if (bunch? (car pls)) (math:error 'elim-bunch? (car pls)))
  27.        (set! npl (cons (poly:elim2 poly (car pls) var) npl))))
  28.       (if (bunch? (car pl)) (math:error 'elim-bunch? (car pl))))))
  29.  
  30. ;(define (poly:restest p1 p2 var)
  31. ;  (let ((v1 (poly:resultant p1 p2 var))
  32. ;    (v2 (poly:elim2 p1 p2 var)))
  33. ;    (cond ((not (equal? v1 v2))
  34. ;       (display-diag "  restest:") (newline-diag)
  35. ;       (math:write (poleqns->licits p1) *output-grammar*)
  36. ;       (newline)
  37. ;       (math:write (poleqns->licits p2) *output-grammar*)
  38. ;       (display-diag "  ==>") (newline-diag)
  39. ;       (math:write (poleqns->licits v1) *output-grammar*)
  40. ;       (display-diag "different from:") (newline-diag)
  41. ;       (math:write (poleqns->licits v2) *output-grammar*)))
  42. ;    v2))
  43.  
  44. (define (intersection? l1 l2)
  45.   (cond ((null? l1) #f)
  46.     ((null? l2) #f)
  47.     ((memq (car l1) l2) #t)
  48.     (else (intersection? (cdr l1) l2))))
  49.  
  50. ;; (sort! ?? (lambda (x y) (< (ext:depth x) (ext:depth y))))
  51. ;(define (ext:depth v)
  52. ;  (let ((vds (var:depends v)))
  53. ;    (if vds (+ 1 (apply max (map ext:depth vds)))
  54. ;    1)))
  55.  
  56. ;;;EVS are all the extension vars used in extensions which are
  57. ;;; not being eliminated.
  58. ;;;IEVS are those EVS which involve VARS.
  59. ;;; sort them so nested extensions are done last.
  60. (define (ext:elim poleqns vars)
  61.   (let* ((eqs (remove-if impl? poleqns))
  62.      (exps (remove-if-not impl? poleqns)))
  63.     (if (> (length exps) 1)
  64.     (math:error 'eliminating-from-more-than-one-expression? exps))
  65.     (let* ((aes (map chainables poleqns))
  66.        (evs (set-difference (reduce union aes) vars))
  67.        (ievs (remove-if-not
  68.           (lambda (ev) (intersection? (var:depends ev) vars))
  69.           evs)))
  70.       (cond
  71.        ((not (null? ievs))
  72.     ;;tievs are the new extensions after any VARS are eliminated
  73.     (do ((ievs ievs (cdr ievs)))
  74.         ((null? ievs))
  75.       (let* ((iev (car ievs))
  76.          (tiev (var:elim iev
  77.                  (remove-if (lambda (x) (poly:find-var? x iev))
  78.                         eqs)
  79.                  vars)))
  80.         (set! eqs (cons (univ:norm0
  81.                  iev
  82.                  (if (expl? tiev) (list tiev -1) (cdr tiev)))
  83.                 eqs))
  84.         (set! vars (cons iev vars))))))
  85.       (poly:elim (append eqs exps) vars))))
  86.  
  87. (define (var:elim var eqs ovars)
  88.   (let* ((vars (var:depends var)))
  89.     (cond ((null? vars) (math:warn 'free-var-to-var:elim var ovars))
  90.       ((radicalvar? var)
  91.        (^ (find-if
  92.            impl?
  93.            (poly:elim (append eqs
  94.                   (list (expl->impl (cadr (func-arglist var)))))
  95.               ovars))
  96.           (caddr (func-arglist var))))
  97.       ((pair? (var:sexp var))
  98.        (apply deferop
  99.           (map (lambda (e)
  100.              (normalize
  101.               (find-if impl? (ext:elim
  102.                       (append eqs (list (licit->impl e)))
  103.                       ovars))))
  104.                (func-arglist var))))
  105.       (else (math:error 'elimination-type-not-handled var)))))
  106.  
  107. (define (infinite-list-of . elts)
  108.   (let ((lst (copy-list elts)))
  109.     (nconc lst lst)))
  110.  
  111. ;;; This tries to solve the equations no matter what is involved.
  112. ;;; It will eliminate variables in vectors of equations.
  113. (define (eliminate eqns vars)
  114.   (bunch:norm
  115.    (if (some bunch? eqns)
  116.        (apply map
  117.           (lambda arglist (eliminate arglist vars))
  118.           (map (lambda (x)
  119.              (if (bunch? x) x (infinite-list-of x)))
  120.            eqns))
  121.        (ext:elim eqns vars))))
  122.  
  123. (define (elim:test)
  124.   (define a (sexp->var 'A))
  125.   (define x (sexp->var 'X))
  126.   (define y (sexp->var 'Y))
  127.   (test (list (list a 0 0 124 81 11 3 45))
  128.     poly:elim
  129.     (list (list y (list x (list a 0 0 2) (list a 0 1)) 1)
  130.           (list y (list x (list a 5 1) 0 -1) 0 1)
  131.           (list y (list x (list a -1 3) 5) -1))
  132.     (list x y)))
  133.  
  134. (define (bunch:map proc b)
  135.   (cond ((bunch? b) (map (lambda (x) (bunch:map proc x)) b))
  136.     (else (proc b))))
  137. (define (licits:for-each proc b)
  138.   (cond ((bunch? b) (for-each (lambda (x) (licits:for-each proc x)) b))
  139.     ((eqn? b) (proc (eqn->poly b)))
  140.     (else (proc b))))
  141. (define (licits:map proc b)
  142.   (cond ((bunch? b) (map (lambda (x) (licits:map proc x)) b))
  143.     ((eqn? b) (poleqn->licit (proc (eqn->poly b))))
  144.     (else (proc b))))
  145. (define (implicits:map proc b)
  146.   (cond ((bunch? b) (map (lambda (x) (implicits:map proc x)) b))
  147.     ((eqn? b) (poleqn->licit (proc (eqn->poly b))))
  148.     ((expl? b) (proc (expl->impl b)))
  149.     (else (proc b))))
  150.  
  151. ;;; replaces each var in poly with (proc var).
  152. ;;; Used for substitutions in clambda and capply.
  153.  
  154. (define (poly:do-vars proc poly)
  155.   (if (number? poly) poly
  156.       (univ:demote (cons (proc (car poly))
  157.              (map (lambda (b) (poly:do-vars proc b))
  158.                   (cdr poly))))))
  159. (define (licits:do-vars proc licit)
  160.   (licits:map (lambda (poly) (poly:do-vars proc poly))
  161.           licit))
  162.  
  163. ;;;; Canonical Lambda
  164. ;;;; This needs to handle algebraic extensions as well.
  165. (define (clambda symlist body)
  166.   (let ((num-new-vars (length (remove-if lambdavar? symlist))))
  167.     (licits:do-vars
  168.      (lambda (var)
  169.        (let ((pos (position (var:nodiffs var) symlist)))
  170.      (cond (pos (lambda-var (+ 1 pos) (var:diff-depth var)))
  171.            ((lambdavar? var) (var:lambda-bump var num-new-vars))
  172.            ((lambdavarext? var) (bump-lambda-ext))
  173.            (else var))))
  174.      body)))
  175.  
  176. (define (clambda? cexp)
  177.   (cond ((number? cexp) #f)
  178.     ((matrix? cexp) (some (lambda (row) (some clambda? row)) cexp))
  179.     ((expr? cexp) (poly:find-var-if? cexp lambdavardep?))
  180.     ((eqn? cexp) (poly:find-var-if? (eqn->poly cexp) lambdavardep?))
  181.     (else #f)))
  182.  
  183. ;;;In order to keep the lambda application hygenic (in case a function
  184. ;;;of a function is called), we need to substitute occurences of
  185. ;;;lambda variables in the body with shadowed versions of the
  186. ;;;variables before we eliminate them.  See:
  187. ;;;    Technical Report No. 194
  188. ;;;    Hygenic Macro Expansion
  189. ;;;    E.E.Kohlbecker, D.P.Friedman, M.Fellinson, and B.Duba
  190. ;;;    Indiana University
  191. ;;;    May, 1986
  192.  
  193. ;;; The bumped-only case is different from the some-bumped
  194. ;;; some-shadowed case in that it returns a publicly available (not
  195. ;;; shadowed) var.  This is called from var:shadow in "types.scm".
  196. (define (var:lambda-bump var delta)
  197.   (if (simple-lambdavar? var)
  198.       (lambda-var (+ (lambda-position var) delta) (var:diff-depth var))
  199.       (sexp->var
  200.        (do-sexp-symbols
  201.     (lambda (s) 
  202.       (define st (symbol->string s))
  203.       (if (and (char=? #\@ (string-ref st 0)) (> (string-length st) 1))
  204.           (var:sexp (lambda-var
  205.              (+ delta (string->number
  206.                    (substring st 1 (string-length st))))
  207.              0))
  208.           s))
  209.     (var:sexp var)))))
  210. (define (do-sexp-symbols proc sexp)
  211.   (cond ((symbol? sexp) (proc sexp))
  212.     ((pair? sexp) (map (lambda (s) (do-sexp-symbols proc s)) sexp))
  213.     (else sexp)))
  214.  
  215. ;;; currently capply puts the structure of the clambda inside the
  216. ;;; structure of the arguments.
  217. (define (capply body larglist)
  218.   (let* ((arglist (licits->impls larglist))
  219.      (arglist-length (length arglist))
  220.      (svlist '())
  221.      (sbody
  222.       (licits:do-vars
  223.        (lambda (var)
  224.          (cond
  225.           ((lambdavardep? var)
  226.            (set! var (var:shadow var arglist-length))
  227.            (set! svlist (union (remove-if-not
  228.                     simple-shadowed-lambdavar?
  229.                     (cons var (var:depends var)))
  230.                    svlist))
  231.            var)
  232.           (else var)))
  233.        body))
  234.      (dargs (diffargs svlist arglist)))
  235.     (implicits:map (lambda (p) (eliminate (cons p dargs) svlist)) sbody)))
  236.  
  237. (define (diffargs vlist args)
  238.   (map (lambda (var)
  239.      (bunch:map (lambda (e)
  240.               (univ:demote (cons var (cdr (licit->impl e)))))
  241.        (diffarg var args)))
  242.     vlist))
  243. (define (diffarg var args)
  244.   (cond ((var:differential? var)
  245.      (total-differential (diffarg (var:undiff var) args)))
  246.     (else (list-ref args (- (lambda-position var) 1)))))
  247.  
  248. (define (licits:for-each-var proc polys)
  249.   (licits:for-each (lambda (poly) (poly:for-each-var proc poly)) polys))
  250.  
  251. (define (licits:max-lambda-position polys)
  252.   (let ((maxpos 0) (deps '()))
  253.     (licits:for-each-var
  254.      (lambda (v) (cond ((lambdavardep? v)
  255.             (set! maxpos (max maxpos (var:max-lambda-position v)))
  256.             (set! deps (adjoin v deps)))))
  257.      polys)
  258.     (for-each
  259.      (lambda (v)
  260.        (for-each
  261.     (lambda (x) (if (lambdavardep? x)
  262.             (set! maxpos (max maxpos (var:max-lambda-position x)))))
  263.     (var:depends v)))
  264.      deps)
  265.     maxpos))
  266.  
  267. (define (var:max-lambda-position var)
  268.   (let ((maxpos 0))
  269.     (for-each
  270.      (lambda (x) 
  271.        (if (lambdavar? x)
  272.        (set! maxpos (max maxpos (or (lambda-position x) maxpos)))))
  273.      (cons var (var:depends var)))
  274.     maxpos))
  275.  
  276. (define (var:min-lambda-position var)
  277.   (let ((minpos 9999))
  278.     (for-each
  279.      (lambda (x) 
  280.        (if (lambdavar? x)
  281.        (set! minpos (min minpos (or (lambda-position x) minpos)))))
  282.      (cons var (var:depends var)))
  283.     (if (= minpos 9999) (error 'no-var:min-lambda-position var))
  284.     minpos))
  285.