home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
ARM Club 3
/
TheARMClub_PDCD3.iso
/
hensa
/
maths
/
b116_1
/
jacal
/
builtin
< prev
next >
Wrap
Text File
|
1993-11-03
|
18KB
|
675 lines
;;; JACAL: Symbolic Mathematics System. -*-scheme-*-
;;; Copyright 1989, 1990, 1991, 1992, 1993 Aubrey Jaffer.
;;; See the file "COPYING" for terms applying to this program.
;;;; First, what case are symbols in? Determine the standard case:
(define char-standard-case
(cond ((not (string=? (symbol->string 'a) (symbol->string 'A)))
char-downcase)
((string=? (symbol->string 'a) "A")
char-upcase)
((string=? (symbol->string 'A) "a")
char-downcase)
(else
char-downcase)))
(define (string-standard-case s)
(set! s (string-copy s))
(do ((i 0 (+ 1 i))
(sl (string-length s)))
((>= i sl) s)
(string-set! s i (char-standard-case (string-ref s i)))))
(define (bltn:error . args)
(apply math:warn args)
novalue)
;;;Predefined Constants
(define expl:t (var->expl (sexp->var 't)))
(define $ (string->var ":@"))
(define $-pri (+ -1 char-code-limit))
(var:set-pri! $ $-pri)
(define ($? v) (or (eq? v $) (= (var:pri v) $-pri)))
(define d$
(var:differential $)) ;used only in total-differential in norm.scm
(var:set-pri! d$ (+ -2 char-code-limit))
(define $1 (string->var "@1"))
(define $2 (string->var "@2"))
(define $3 (string->var "@3"))
(define _$ (string->var "::@"))
(define $1+$2 (list $2 (list $1 0 1) 1))
(define $1-$2 (list $2 (list $1 0 1) -1))
(define _-$1 (list $1 0 -1))
(define $1*$2 (list $2 0 (list $1 0 1)))
(define $1/$2 (list $ (list $1 0 1) (list $2 0 -1)))
(define $1=$2 (list $= $2 (list $1 0 -1) 1))
(define cidentity (list $1 0 1))
(define _^ (list (string->var "^") 0 1))
(define _^^ (list (string->var "^^") 0 1))
(define _partial (list (string->var "partial") 0 1))
(define _ncmult (list (string->var "ncmult") 0 1))
;;; canoncial functions for vect.scm
(define $1-$2*$3 (list $3 (list $1 0 1) (list $2 0 -1)))
(define _-$1/$2 (make-rat (list $1 0 -1) (list $2 0 1)))
(define $1*$2+$3 (list $3 (list $2 0 (list $1 0 1)) 1))
;;; set up initial radical and extension
(define %sqrt1 (defext (sexp->var '%sqrt1) (list $ 1 0 -1)))
(var:set-pri! %sqrt1 5)
(define %i (defext (sexp->var '%i) (list $ -1 0 -1)))
(var:set-pri! %i 5)
(define radical-defs (list (extrule %i) (extrule %sqrt1)))
(define _+/-$1 (list $1 0 (list %sqrt1 0 1)))
(define _-/+$1 (list $1 0 (list %sqrt1 0 -1)))
(define $1+/-$2 (list $2 (list $1 0 1) (list %sqrt1 0 1)))
(define $1-/+$2 (list $2 (list $1 0 1) (list %sqrt1 0 -1)))
(define novalue (var->expl (sexp->var '?)))
(define (novalue? x) (equal? novalue x))
(define *flags* '())
(define flag-associator (alist-associator eq?))
(define flag-inquirer (alist-inquirer eq?))
(define (list-of-flags)
(define flags '())
(alist-for-each (lambda (k v) (set! flags (cons k flags))) *flags*)
flags)
;(define *flags* (make-hash-table 5))
;(define flag-associator (hash-associator eq?))
;(define flag-inquirer (hash-inquirer eq?))
;(define (list-of-flags)
; (define flags '())
; (hash-for-each (lambda (k v) (set! flags (cons k flags))) *flags*)
; flags)
(define (defflag name setter getter)
(set! *flags* (flag-associator *flags* name (cons setter getter)))
name)
(define flag:setter car)
(define flag:getter cdr)
(define (flag-set name . values)
(let ((flag (flag-inquirer *flags* name)))
(cond ((not flag) (bltn:error 'flag name 'is-not-defined))
((flag:setter flag) (apply (flag:setter flag) flag values) novalue)
(else (bltn:error 'flag name 'can-not-be-set)))))
(define (flag-get name . rest)
(let ((flag (flag-inquirer *flags* name)))
(cond ((not flag) (bltn:error 'flag name 'is-not-defined))
((flag:getter flag) (apply (flag:getter flag) flag rest))
(else (bltn:error 'flag name 'can-not-be-read)))))
(defflag 'ingrammar
(lambda (f v)
(define name (var:sexp (expl->var v)))
(cond ((get-grammar name)
(set! *input-grammar* (get-grammar name)))
(else
(bltn:error 'grammar name 'not-known))))
(lambda (f) (var->expl (sexp->var (grammar-name *input-grammar*)))))
(defflag 'outgrammar
(lambda (f v)
(define name (var:sexp (expl->var v)))
(cond ((get-grammar name)
(set! *output-grammar* (get-grammar name)))
(else
(bltn:error 'grammar name 'not-known))))
(lambda (f) (var->expl (sexp->var (grammar-name *output-grammar*)))))
(defflag 'echogrammar
(lambda (f v)
(define name (var:sexp (expl->var v)))
(cond ((get-grammar name)
(set! *echo-grammar* (get-grammar name)))
(else
(bltn:error 'grammar name 'not-known))))
(lambda (f) (var->expl (sexp->var (grammar-name *echo-grammar*)))))
(defflag 'grammars
#f
(lambda (f)
(map (lambda (g) (var->expl (sexp->var g))) (list-of-grammars))))
(define (set-boolean v)
(define val (var:sexp (expl->var v)))
(case val
((off 0 false) #f)
((on 1 true) #t)
(else (bltn:error 'expected-boolean v))))
(define (show-boolean v)
(var->expl (sexp->var (if v 'on 'off))))
(defflag 'horner
(lambda (f v) (set! horner (set-boolean v)))
(lambda (f) (show-boolean horner)))
(defflag 'trace
(lambda (f v) (set! math:trace (set-boolean v)))
(lambda (f) (show-boolean math:trace)))
(defflag 'debug
(lambda (f v) (set! math:debug (set-boolean v)))
(lambda (f) (show-boolean math:debug)))
(defflag 'phases
(lambda (f v) (set! math:phases (set-boolean v)))
(lambda (f) (show-boolean math:phases)))
(defflag 'linkradicals
(lambda (f v) (set! linkradicals (set-boolean v)))
(lambda (f) (show-boolean linkradicals)))
(defflag 'version
#f
(lambda (f)
(var->expl (string->var *jacal-version*))))
(defflag 'all
#f
(lambda (f)
(block-write-strings
(sort! (map symbol->string (list-of-flags))
string<?))
novalue))
(defflag 'prompt
(lambda (f v)
(set! newlabelstr (var->string (expl->var v)))
(set! newlabelsym (string->symbol newlabelstr))
novalue)
(lambda (f) (var->expl (string->var newlabelstr))))
(defflag 'page
(lambda (f v)
(define val (if (number? v) v (var:sexp (expl->var v))))
(set! page-height
(case val ((off 0 false) #f)
((on 1 true) #t)
(else (if (number? val) val
(bltn:error 'expected-boolean-or-number v))))))
(lambda (f) (if (boolean? page-height)
(show-boolean page-height)
page-height)))
(defflag 'width
(lambda (f v)
(define val (if (number? v) v (var:sexp (expl->var v))))
(set! page-width
(case val ((off 0 false) #f)
((on 1 true) #t)
(else (if (number? val) val
(bltn:error 'expected-boolean-or-number v))))))
(lambda (f) (if (boolean? page-width)
(show-boolean page-width)
page-width)))
(defflag 'priority
(lambda (f v p)
(if (not (and (number? p) (< 0 p lambda-var-pri))) (math:error))
(var:set-pri! (expl->var v) p))
(lambda args
(if (null? (cdr args))
(let ((l (list-of-vars)))
(block-write-strings (map object->string
(map var:sexp (sort! l var:>))))
novalue)
(var:pri (expl->var (cadr args))))))
;(define transcript-name #f)
;(defflag 'transcript
; (lambda (f v)
; (define file (and v (not (null? v)) (var->string (expl->var v))))
; (if v (transcript-on file) (transcript-off))
; (set! transcript-name file))
; (lambda (f) (if transcript-name
; (var->expl (string->var transcript-name))
; '#())))
;;;; Built in functions
(defbltn 'set
(lambda (name . values)
(apply flag-set (var:sexp (expl->var name)) values)))
(defbltn 'show
(lambda (name . rest) (apply flag-get
(var:sexp (expl->var name))
rest)))
(defbltn 'commands
(lambda ()
(block-write-strings
(sort! (map object->string (list-of-procedure-defsyms))
string<?))
novalue))
(defbltn '%
(lambda () %))
(defbltn 'depends
(lambda (x) (map var->expl (var:depends (expl->var x)))))
(defbltn 'args
(lambda (x) (cdr (func-arglist (expl->var x)))))
(defbltn 'func
(lambda (x) (car (func-arglist (expl->var x)))))
(defbltn 'describe
(lambda (x)
(cond
((and (expl:var? x)
(info:describe (var:sexp (expl->var x)))))
((bunch? x) (display (bunch-type x)) (newline))
((not (expl:var? x)) (display (scalar-type x)))
(else (describe-var (expl->var x))))
(if (clambda? x)
(let ((hlv (licits:max-lambda-position (if (eqn? x) (eqn->poly x) x))))
(tran:display 'function-of-)
(display hlv)
(if (= 1 hlv) (tran:display 'argument) (tran:display 'arguments))))
novalue))
(define (describe-var v)
(cond ((var:differential? v)
(tran:display 'differential-)
(set! v (var:nodiffs v))))
(display
(cond ((radicalvar? v) 'radical)
((not (symbol? (var:sexp v))) 'application)
(else 'variable))))
(define (scalar-type x)
(cond ((number? x) 'number)
((eqn? x) 'equation)
((expl? x) 'polynomial)
((rat? x) 'rational-expression)
((impl? x) 'implicit-expression)
(else 'unknown)))
(define (bunch-type x)
(cond ((matrix? x) 'matrix)
((row? x) 'row-vector)
((column? x) 'column-vector)
(else 'bunch)))
(defbltn 'example
(lambda (x) (info:example x)))
(define (terms) (paginate-file (in-vicinity jacal-vicinity "COPYING")))
(defbltn 'terms (lambda () (terms) novalue))
(define (help) (paginate-file (in-vicinity jacal-vicinity "HELP")))
(defbltn 'help (lambda () (help) novalue))
(define (boolify x)
(var->expl (sexp->var (if x 'true 'false))))
(defbltn 'verify
(lambda (try expect)
(let ((tv (normalize try)) (ev (normalize expect)))
(cond ((equal? tv ev) (boolify #t))
(else
(display-diag (tran:translate 'Did-not-verify:))
(newline-diag)
(write-sexp (math->sexp (list tv ev) horner) *output-grammar*)
(if math:debug (do-more))
(boolify #f))))))
(defbltn 'Differential
(lambda (obj) (total-differential obj)))
(defbltn 'negate
(lambda (obj) (app* _-$1 obj)))
(defbltn 'u+/-
(lambda (obj) (app* _+/-$1 obj)))
(defbltn 'u-/+
(lambda (obj) (app* _-/+$1 obj)))
(defbltn '^ ;need to do expt also
(lambda (x exp)
(if (and (expl? x) (number? exp) (positive? exp))
(poly:^ x (normalize exp))
(^ (expr x) exp))))
(defbltn '^^ ;need to do ncexpt also
(lambda (a pow) (ncexpt (exprs a) (normalize pow))))
(defbltn '*
(lambda args (reduce (lambda (x y)
(if (and (expl? x) (expl? y))
(poly:* x y)
(app* $1*$2 x y)))
args)))
(defbltn '+
(lambda args (reduce (lambda (x y)
(if (and (expl? x) (expl? y))
(poly:+ x y)
(app* $1+$2 x y)))
args)))
(defbltn '-
(lambda args (reduce (lambda (x y)
(if (and (expl? x) (expl? y))
(poly:- x y)
(app* $1-$2 x y)))
args)))
(defbltn 'b+/-
(lambda args (reduce (lambda (x y) (app* $1+/-$2 x y)) args)))
(defbltn 'b-/+
(lambda args (reduce (lambda (x y) (app* $1-/+$2 x y)) args)))
(defbltn '/
(lambda args (reduce (lambda (x y) (app* $1/$2 x y)) args)))
(defbltn 'over
(lambda args (reduce (lambda (x y) (app* $1/$2 x y)) args)))
(defbltn 'bunch
(lambda args args))
(defbltn 'rapply
(lambda args (apply rapply args)))
(defbltn 'or
(lambda args
(poleqn->licit (reduce poly:* (map licit->poleqn args)))))
(defbltn '=
(lambda (x y) (app* $1=$2 x y)))
(defbltn 'qed
(lambda ()
(cleanup-handlers!)
(math:exit #t)))
(defbltn 'quit
(lambda ()
(cleanup-handlers!)
(quit)))
;;;; User callable functions
(defbltn 'listofvars
(lambda (exp)
(let ((deps '()))
(licits:for-each (lambda (poly) (set! deps (union (alg:vars poly) deps)))
exp)
(map var->expl (remove $ deps)))))
(defbltn 'coeff
(lambda (p var . optional)
(let ((ord (if (null? optional) 1 (car optional))))
(poly:coeff p (expl->var var) (plicit->integer ord)))))
(defbltn 'num
(lambda (exp) (num (expr:normalize exp))))
(defbltn 'denom
(lambda (exp) (denom (expr:normalize exp))))
(defbltn 'divide
(lambda (dividend divisor . vars)
(set! dividend (licit->polxpr dividend))
(set! divisor (licit->polxpr divisor))
(poly:pdiv dividend divisor (if (null? vars)
(if (number? divisor)
(if (number? dividend) 0
(car dividend))
(car divisor))
(expl->var (car vars))))))
(defbltn 'content
(lambda (poly var)
(let* ((var (expl->var var))
(poly (promote var (licit->polxpr poly)))
(cont (apply poly:gcd* (cdr poly))))
(list cont (poly:/ poly cont)))))
;;; This is user callable GCD.
(defbltn 'gcd
(lambda args
(if (null? args) 0
(reduce poly:gcd (map licit->polxpr args)))))
(defbltn 'mod
(lambda (licit polxpr)
(poleqn->licit (alg:mod (licit->poleqn licit) (licit->polxpr polxpr)))))
;;; This is user callable RESULTANT. It always operates on
;;; polynomials and does not know about extensions etc.
(defbltn 'resultant
(lambda (a b v)
(let ((res (poly:resultant
(licit->polxpr a)
(licit->polxpr b)
(expl->var v))))
res)))
(defbltn 'sylvester
(lambda (p1 p2 var)
(sylvester (licit->polxpr p1)
(licit->polxpr p2)
(expl->var var))))
(defbltn 'discriminant
(lambda (poly var)
(set! poly (licit->polxpr poly))
(set! poly (poly:/ poly (if (> (leading-number poly) 0)
(poly:num-cont poly)
(- (poly:num-cont poly)))))
(let* ((v (expl->var var))
(deg (poly:degree poly v)))
(if (> deg 1)
(poly:* (quotient (* deg (- deg 1)) 2)
(poly:resultant (poly:diff poly v) poly v))
0))))
(defbltn 'eliminate
(lambda (eqns vars)
(poleqns->licits (eliminate (licits->poleqns eqns) (variables (normalize vars))))))
(defbltn 'PolyElim
(lambda (eqns vars)
(poleqns->licits (poly:elim (licits->poleqns eqns) (variables vars)))))
(defbltn 'factor
(lambda (poly)
(let ((e (licit->polxpr poly)))
(cond ((number? e) (require 'prime) ;autoload from SLIB
(sort! (factor e) <))
(else (poly:factorq e))))))
(defbltn 'prime?
(lambda (n)
(let ((e (licit->polxpr n)))
(cond ((number? e) (require 'prime) ;autoload from SLIB
(boolify (prime? e)))
(else (bltn:error 'Not-a-Number n))))))
(defbltn 'matrix
(lambda args (apply matrix args)))
(defbltn 'genmatrix
(lambda (fun i2 j2 . i1j1)
(let ((i1 1) (j1 1))
(cond ((null? i1j1))
((begin (set! i1 (car i1j1))
(set! i1j1 (cdr i1j1))
(set! j1 i1)
(null? i1j1)))
((begin (set! j1 (car i1j1))
(set! i1j1 (cdr i1j1))
(null? i1j1)))
(else (math:error 'genmatrix wna)))
(mtrx:genmatrix
fun
(plicit->integer i2)
(plicit->integer j2)
(plicit->integer i1)
(plicit->integer j1)))))
(defbltn 'augcoefmatrix
(lambda (eqns vars)
(augcoefmatrix (licits->poleqns eqns) (variables vars))))
(defbltn 'coefmatrix
(lambda (eqns vars)
(coefmatrix (licits->poleqns eqns) (variables vars))))
(defbltn 'rank
rank)
(defbltn 'ident
(lambda (n) (mtrx:scalarmatrix n 1)))
(defbltn 'scalarmatrix
(lambda (n x) (mtrx:scalarmatrix (plicit->integer n) x)))
(defbltn 'diagmatrix
(lambda args (mtrx:diagmatrix args)))
(defbltn 'determinant
(lambda (m) (determinant m)))
(defbltn 'charpoly
charpoly)
(defbltn 'crossproduct
(lambda (x y) (crossproduct x y)))
(defbltn 'dotproduct
(lambda (x y) (dotproduct x y)))
(defbltn 'ncmult
(lambda (x y) (ncmult x y)))
(defbltn 'row
(lambda (m i)
(if (matrix? m)
(list-ref m (+ -1 (plicit->integer i)))
(bltn:error 'Row-of-non-matrix?:- M))))
(defbltn 'col
(lambda (m i)
(cond ((matrix? m)
(map (lambda (row)
(list (list-ref row (+ -1 (plicit->integer i)))))
m))
((bunch? m) (list-ref m (plicit->integer i)))
(else (bltn:error 'Column-of-non-matrix?:- M)))))
(defbltn 'minor
(lambda (m i j)
(mtrx:minor m (plicit->integer i) (plicit->integer j))))
(defbltn 'cofactor
(lambda (m i j)
(cofactor m (plicit->integer i) (plicit->integer j))))
(defbltn 'transpose
(lambda (m) (transpose m)))
(defbltn 'CartProd
(lambda (m) (cart-prod m)))
(defbltn 'Elementwise
(lambda (f . args)
(apply map (lambda args (sapply f args)) args)))
(defbltn 'finv
(lambda (f)
(fcinverse f)))
(defbltn 'load
(lambda (file)
(load (var->string (expl->var file)))
file))
(defbltn 'require
(lambda (file)
(load (in-vicinity jacal-vicinity (var->string (expl->var file))))
file))
(defbltn 'batch
(lambda (file)
(batch (var->string (expl->var file)))
novalue))
(defbltn 'transcript
(lambda files
(cond ((null? files)
(transcript-off)
novalue)
((not (null? (cdr files))) (bltn:error 'transcript wna files))
(else
(let ((file (var->string (expl->var (car files)))))
(transcript-on file)
(car files))))))
(defbltn 'system
(lambda (command)
(system (var->string (expl->var command)))
; command ;uncomment this line if system doesn't return nicely
))
(defbltn 'coeffs
(lambda (poly var)
(if (not (and (expl? poly) (not (number? poly))))
(bltn:error 'not-a-polynomial? poly)
(cdr (promote (expl->var var) poly)))))
(defbltn 'poly
(lambda (var . args)
(reduce (lambda (p c) (poly:+ (poly:* p var) c))
(cond ((> (length args) 1) args)
((not (and (= (length args) 1) (bunch? (car args))))
(bltn:error 'not-a-bunch? (car args)))
(else (car args))))))
(defbltn 'diff
(lambda (exp . args)
(reduce-init diff exp (map expl->var args))))
(defbltn 'PolyDiff
(lambda (exp . args)
(reduce-init expls:diff exp (map expl->var args))))
(defbltn 'partial
(lambda (func . args)
(cond ((number? func) (bltn:error 'not-a-function? func))
((null? args) (bltn:error 'no-variables?))
((not (clambda? func)) (apply deferop _partial func args))
(else
(reduce-init
diff func
(map (lambda (a)
(cond ((and (number? a) (positive? a)) (lambda-var a 0))
((clambda? a) (expl->var a))
(else (math:error 'partial-with-respect-to? a))))
args))))))
;;; commands for debugging:
(defbltn 'chain
(lambda (exp)
(let ((e (expl->var exp)))
(poly->eqn (chain-rule e (var:differential e))))))
(defbltn 'shadow
(lambda (x) (map (lambda (v) (if v (var->expl v) '()))
(or (vector-ref (expl->var x) 4) '()))))
(defbltn 'extrule
(lambda (x) (poly->eqn (or (extrule (expl->var x)) 0))))