;;
;; Definition and code generators for a simple applicative language
;; for numerical simulation.
;;
;; Copyright 2010-2013 Ivan Raikov and the Okinawa Institute of
;; Science and Technology.
;;
;; This program is free software: you can redistribute it and/or
;; modify it under the terms of the GNU General Public License as
;; published by the Free Software Foundation, either version 3 of the
;; License, or (at your option) any later version.
;;
;; This program is distributed in the hope that it will be useful, but
;; WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
;; General Public License for more details.
;;
;; A full copy of the GPL license can be found at
;; .
;;
(module flsim
(
value?
V:C V:Var V:Rec V:Sel V:Vec V:Sub V:Fn V:Op V:Ifv
expr?
E:Ife E:Let E:Ret E:Noop
binding? B:Val
typenv
name/scheme prelude/scheme expr->scheme value->scheme binding->scheme
name/ML prelude/ML expr->ML value->ML binding->ML
name/Octave prelude/Octave expr->Octave value->Octave binding->Octave
)
(import scheme chicken)
(require-extension extras data-structures srfi-1 datatype)
(define nl "\n")
(define-datatype value value?
(V:C (v (lambda (x) (or (symbol? x) (number? x) ))))
(V:Var (name symbol?))
(V:Rec (flds (lambda (xs) (and (pair? xs) (every (lambda (x) (and (symbol? (car x)) (value? (cadr x)))) xs)))))
(V:Sel (field (lambda (x) (or (symbol? x) (integer? x)))) (v value?))
(V:Vec (vals (lambda (x) (every value? x))))
(V:Sub (index (lambda (x) (and (integer? x) (or (zero? x) (positive? x))))) (v value?))
(V:Fn (args list?) (body expr?))
(V:Op (name symbol?) (args (lambda (x) (every value? x))))
(V:Ifv (test value?) (ift value?) (iff value?))
)
(define-record-printer (value x out)
(fprintf out "~A"
(cases value x
(V:C (v) (sprintf "V:C ~A" v))
(V:Var (n) (sprintf "V:Var ~A " n))
(V:Rec (lst) (sprintf "V:Rec ~A" lst))
(V:Sel (f v) (sprintf "V:Sel ~A ~A" f v))
(V:Vec (lst) (sprintf "V:Vec ~A" lst))
(V:Sub (i v) (sprintf "V:Sub ~A ~A" i v))
(V:Ifv (test tv fv) "V:Ifv ~A ~A ~A" test tv fv)
(V:Fn (args body) (sprintf "V:Fn ~A = ~A" args body))
(V:Op (name args) (sprintf "(V:Op ~A ~A)" name args)))))
(define-datatype binding binding?
(B:Val (name symbol?) (v value?))
)
(define-record-printer (binding x out)
(fprintf out "~A"
(cases binding x
(B:Val (name v) (sprintf "B:Val ~A = ~A" name v))
)))
(define-datatype expr expr?
(E:Ife (test value?) (ift expr?) (iff expr?))
(E:Let (bnds (lambda (x) (every binding? x))) (body expr?))
(E:Ret (v value?))
(E:Noop)
)
(define-record-printer (expr x out)
(fprintf out "~A"
(cases expr x
(E:Ife (test ift iff) (sprintf "E:Ife ~A ~A ~A" test ift iff))
(E:Let (bnds body) (sprintf "E:Let ( ~A ) ~A" bnds body))
(E:Ret (v) (sprintf "E:Ret ~A" v))
(E:Noop () (sprintf "E:Noop"))
)))
;; Builds a type environment for an expression list
(define (typenv bindings)
(define (numeric-type? x)
(case x ((numeric) #t) (else #f)))
(define (boolean-type? x)
(case x ((bool) #t) (else #f)))
(define (type-equal? ty1 ty2)
(case (car ty1)
((numeric boolean unit)
(eq? (car ty1) (car ty2)))
((record)
(case (car ty2)
((record)
(lset= (lambda (x y) (and (eq? (car x) (car y))
(type-equal? (cdr x) (cdr y))))
(cdr ty1) (cdr ty2)))
(else #f)))
((vector)
(case (car ty2)
((vector)
(and (= (length (cdr ty1)) (length (cdr ty2)))
(every type-equal? (cdr ty1) (cdr ty2))))
(else #f)))
((function)
(case (car ty2)
((function)
(and (equal? (cadr ty1) (cadr ty2))
(equal? (caddr ty1) (caddr ty2))))))
))
(define (value-type v env)
(cases value v
(V:C (v) (if (number? v) 'number
(case v
((false true) 'bool)
(else (error 'value-type "unknown constant" v)))))
(V:Var (s) (alist-ref s env))
(V:Rec (flds)
`(record . ,(map (lambda (f) (cons (car f) (value-type (cadr f) env))) flds)))
(V:Sel (field u)
(let ((ut (value-type u env)))
(cond ((and (eq? (car ut) 'record) (symbol? field))
(alist-ref field (cdr ut)))
(else (error 'value-type "invalid select value" field u))
)))
(V:Vec (vs) `(vector . ,(map (lambda (x) (value-type x env)) vs)))
(V:Sub (index u)
(let ((ut (value-type u env)))
(cond
((and (eq? (car ut) 'vector) (integer? index))
(let ((limit (length (cdr ut))))
(if (not (or (zero? index) (positive? index)))
(error 'value-type "invalid vector index" index u))
(if (< index limit)
(error 'value-type "vector index out of range" index u))
(list-ref (cdr ut) index)
))
(else (error 'value-type "invalid sub value" index u))
)))
(V:Fn (args body) `(function ,args ,body))
(V:Op (name args)
(let ((f (alist-ref name env)))
(if f
(apply-type f args env)
(case name
((+ - * /)
(if (every numeric-type? (map (lambda (x) (value-type x env)) args))
'number
(error 'value-type "non-numeric arguments to arithmetic operator"
name args)))
((>= > < <=)
(if (every numeric-type? (map (lambda (x) (value-type x env)) args))
'bool
(error 'value-type "non-numeric arguments to comparison operator"
name args)))
((TRSA TRSB)
(let ((ty (value-type (car args) env)))
`(record (,name . ,ty))))
(else (error 'value-type "unknown operator" name))
))
))
(V:Ifv (test ift iff)
(if (boolean-type? (value-type test env))
(let ((ty1 (value-type ift env))
(ty2 (value-type iff env)))
(if (type-equal? ty1 ty2) ty1
(error 'value-type "if value branches have different types"
test ift ty1 iff ty2)))
(error 'value-type "if value has non-boolean test value")))
))
(define (binding-env b env)
(cases binding b
(B:Val (name v)
(let ((vt (value-type v env)))
(cons `(name . vt) env)))
))
(define (expr-type e env)
(cases expr e
(E:Ife (test ift iff)
(if (boolean-type? (value-type test env))
(let ((ty1 (expr-type ift env))
(ty2 (expr-type iff env)))
(if (type-equal? ty1 ty2) ty1
(error 'expr-type "if expression branches have different types"
test ift ty1 iff ty2)))
(error 'expr-type "if expression has non-boolean test value")))
(E:Let (bnds body)
(let ((env1 (fold binding-env env bnds)))
(expr-type e env1)))
(E:Ret (v) (value-type v env))
(E:Noop () 'unit)))
(define (apply-type f args env)
(let ((formals (cadr f))
(body (caddr f)))
(let ((env (map (lambda (formal arg) `(,formal . ,(value-type arg env))) formals args)))
(expr-type body env))))
(fold binding-env '() bindings)
)
(define (name/scheme s)
(let ((cs (string->list (->string s))))
(let loop ((lst (list)) (cs cs))
(if (null? cs) (string->symbol (list->string (reverse lst)))
(let* ((c (car cs))
(c1 (cond ((or (char-alphabetic? c) (char-numeric? c)
(char=? c #\_) (char=? c #\-)) c)
(else #\-))))
(loop (cons c1 lst) (cdr cs)))))))
(define (binding->scheme x)
(cases binding x
(B:Val (name v)
(list "(" (name/scheme name) " " (value->scheme v) ")" nl))))
(define (expr->scheme x . rest)
(let-optionals rest ((bnd? #t))
(cases expr x
(E:Ife (test ift iff)
(list "(cond (" (value->scheme test) " " nl
" " (expr->scheme ift ) ")" nl
"(else " (expr->scheme iff) "))" nl))
(E:Let (bnds body)
(list "(let* (" nl
(map (lambda (x) (binding->scheme x)) bnds) nl
") " nl
(expr->scheme body #f) nl
")" nl))
(E:Ret (v) (value->scheme v))
(E:Noop () (list "(void)"))
)))
(define (value->scheme v)
(let ((result
(cases value v
(V:C (v) v)
(V:Var (name) (name/scheme name))
(V:Rec (lst)
(list "`(" (intersperse (map (lambda (nv) (list "(" (name/scheme (car nv)) " . ,"
(value->scheme (cadr nv)) ")")) lst) " ") ")"))
(V:Sel (field v)
(if (number? field)
(list "(list-ref " (value->scheme v) " " (- field 1) ")")
(list "(alist-ref '" (name/scheme field) " " (value->scheme v) ")")))
(V:Vec (lst)
(list "(list " (intersperse (map value->scheme lst) " ") ")"))
(V:Sub (index v)
(list "(list-ref " (value->scheme v) " " index ")"))
(V:Fn (args body)
(list "(lambda (" (intersperse (map name/scheme args) " ") ") "
(expr->scheme body #f) ")"))
(V:Op (name args)
(let* ((fp? (case name
((+ - * / >= > < <= neg) #t)
(else #f)))
(op (if fp? (conc "fp" name) name)))
(cond ((null? args)
(case name
((NONE) (list "#f"))
(else (list "(" name ")"))))
(fp?
(if (pair? (cdr args))
(fold-right (lambda (x ax) (list "(" op " " (value->scheme x) " " ax ")"))
(list "(" op " " (value->scheme (car args)) " " (value->scheme (cadr args)) ")")
(cddr args))
(list "(" op " " (value->scheme (car args)) ")")
))
(else
(list "(" op " " (intersperse (map value->scheme args) " ") ")")))))
(V:Ifv (test ift iff)
(list "(if " (value->scheme test) " "
(value->scheme ift) " "
(value->scheme iff) ")"))
)))
result))
(define (prelude/scheme #!key (solver 'rk4b) (integral-index 0))
`(
,(case solver
((cvode) `("(use sundials random-mtzig mathh)" ,nl))
(else `("(use runge-kutta random-mtzig mathh)" ,nl)))
#<list (->string s))))
(let loop ((lst (list)) (cs cs))
(cond ((null? cs) (string->symbol (list->string (reverse lst))))
((null? (cdr cs))
(let ((c (car cs)))
(if (or (char-alphabetic? c) (char-numeric? c))
(loop (cons c lst) (cdr cs))
(loop (append (reverse (string->list (->string (gensym 't)))) lst) (cdr cs))
)))
(else
(let* ((c (car cs))
(c1 (cond ((or (char-alphabetic? c) (char-numeric? c)
(char=? c #\_) (char=? c #\#)) c)
(else #\_))))
(loop (cons c1 lst) (cdr cs))))))))
(define (binding->ML x . rest)
(let-optionals rest ((bnd? #t))
(cases binding x
(B:Val (name v)
(list "val " (name/ML name) " = " (value->ML v) nl))
)))
(define (expr->ML x . rest)
(let-optionals rest ((bnd? #t))
(cases expr x
(E:Ife (test ift iff)
(list "if (" (value->ML test) ") " nl
"then " (expr->ML ift ) nl
"else " (expr->ML iff) nl))
(E:Let (bnds body)
(list "let " nl
(map (lambda (x) (binding->ML x #t)) bnds) nl
"in " nl
(expr->ML body #f) nl
"end" nl))
(E:Ret (v) (value->ML v))
(E:Noop () (list "()"))
)))
(define (value->ML v)
(cases value v
(V:C (v) (cond ((and (number? v) (negative? v))
(list "~" (abs v)))
(else v)))
(V:Var (name) (name/ML name))
(V:Rec (lst)
(list "{" (intersperse (map (lambda (nv) (list (name/ML (first nv)) " = "
(value->ML (cadr nv)))) lst) ", ") "}"))
(V:Sel (field v)
(list "(#" (name/ML field) "(" (value->ML v) "))"))
(V:Vec (lst)
(let ((n (length lst)))
(list "([" (intersperse (map (lambda (v) (value->ML v)) lst) ", ") "])")))
(V:Sub (index v)
(list "List.nth (" (value->ML v) ", " index ")"))
(V:Fn (args body)
(list "(fn (" (intersperse (map name/ML args) ",") ") => "
(expr->ML body #f) ")"))
(V:Op (name args)
(let* ((infix? (case name
((+ - * / >= > < <=) #t)
(else #f)))
(op (if infix? (list "(op " name ")") name)))
(cond ((null? args)
(case name
((NONE) (list name))
(else (list name "()"))))
((null? (cdr args))
(list "(" op " " (value->ML (car args)) ")"))
((and infix? (pair? (cddr args)))
(list "(foldr " op "(" (value->ML (V:Op name (list (car args) (cadr args)))) ")"
"[" (intersperse (map value->ML (cddr args)) ",") "])"))
(else
(list "(" op "(" (intersperse (map value->ML args) ",") "))")))))
(V:Ifv (test ift iff)
(list "(if (" (value->ML test) ") "
"then " (value->ML ift ) " "
"else " (value->ML iff) ")"))
))
(define (prelude/ML #!key (solver 'rk4b) )
`(
#< (('b,'c) trs))) *
(('a -> (('b,'c) trs))) *
((('b,'c) trs) -> bool))
fun tsCase (fa,fb,x) = case x of TRSA a => (fa a) | TRSB b => (fb b)
fun trfOf x = case x of TRC (f,fk,e) => f
fun trfkOf x = case x of TRC (f,fk,e) => fk
fun treOf x = case x of TRC (f,fk,e) => e
fun putStrLn str =
(TextIO.output (TextIO.stdOut, str);
TextIO.output (TextIO.stdOut, "\n"))
fun putStr str = (TextIO.output (TextIO.stdOut, str))
fun showReal n =
let open StringCvt
in
(if n < 0.0 then "-" else "") ^ (fmt (FIX (SOME 12)) (abs n))
end
exception EmptySignal
val neg = (op ~)
val swap = fn (x,v) => (case v of NONE => x | SOME v => v)
val equal = fn (x,y) => (x = y)
val signalOf = fn (v) => (case v of NONE => raise EmptySignal | SOME v => v)
fun PoissonInit () =
let
val zt = RandomMTZig.ztnew()
val st = RandomMTZig.fromEntropy()
in
{st=st,zt=zt}
end
fun PoissonStep (rate,t,h,st,zt) =
let
val rv = RandomMTZig.randPoisson (rate*0.001*h,st,zt)
val spike' = Real.> (rv,0.0)
val spikeCount' = if spike' then rv else 0.0
in
{t=t+h,spike=spike',spikeCount=spikeCount',st=st,zt=zt}
end
EOF
,(if (not solver)
`(("(* dummy solver; returns only the computed derivatives *)")
("fun integrate (f,x: real,y: real list,h,i) = (f (x,y))" ,nl)
)
`(("val summer = fn (a,b) => (ListPair.map (fn (x,y) => x+y) (a,b))" ,nl)
("val scaler = fn(a,lst) => (map (fn (x) => a*x) lst)" ,nl)
("val " ,solver ": (real list) stepper1 = make_" ,solver "()" ,nl)
("fun make_stepper (deriv) = " ,solver " (scaler,summer,deriv)" ,nl)
("fun integrate (f,x: real,y: real list,h,i) = ((make_stepper f) h) (x,y)" ,nl)
))
))
(define (name/Octave s)
(let ((cs (string->list (->string s))))
(let loop ((lst (list)) (cs cs))
(cond ((null? cs) (string->symbol (list->string (reverse lst))))
((null? (cdr cs))
(let ((c (car cs)))
(if (or (char-alphabetic? c) (char-numeric? c))
(loop (cons c lst) (cdr cs))
(loop (append (reverse (string->list (->string (gensym 't)))) lst) (cdr cs))
)))
(else
(let* ((c (car cs))
(c1 (cond ((or (char-alphabetic? c) (char-numeric? c)
(char=? c #\_)) c)
(else #\_))))
(loop (cons c1 lst) (cdr cs))))))))
(define (binding->Octave x . rest)
(let-optionals rest ((inner? #f))
(cases binding x
(B:Val (name v)
(if inner?
(list (name/Octave name) " = " (value->Octave v) )
(list (name/Octave name) " = " (value->Octave v) ";" nl))))))
(define (expr->Octave x . rest)
(let-optionals rest ((inner? #f))
(cases expr x
(E:Ife (test ift iff)
(list "if (" (value->Octave test) ") " nl
(expr->Octave ift ) nl
" else " (expr->Octave iff) " endif" nl))
(E:Let (bnds body)
(if inner?
(list #\{
(intersperse (append (map (lambda (x) (list "(" (binding->Octave x #t) ")")) bnds)
(list "(" (expr->Octave body #t) ")"))
", ")
#\}
"{" (+ 1 (length bnds)) "}" nl)
(list (map (lambda (x) (list (binding->Octave x) nl)) bnds) nl
(expr->Octave body) nl)))
(E:Ret (v) (value->Octave v))
(E:Noop () (list))
)))
(define (value->Octave v)
(cases value v
(V:C (v) v)
(V:Var (name) (name/Octave name))
(V:Rec (lst)
(list "struct (" (intersperse (map (lambda (nv) (list #\" (name/Octave (first nv)) #\" ", "
(value->Octave (cadr nv)))) lst) ", ") ")"))
(V:Sel (field v)
(list (value->Octave v) "." (name/Octave field)))
(V:Vec (lst)
(let ((n (length lst)))
(list "([" (intersperse (map (lambda (v) (value->Octave v)) lst) ", ") "])")))
(V:Sub (index v)
(list (value->Octave v) "((" index ")+1" ")"))
(V:Fn (args body)
(list "(@ (" (intersperse (map name/Octave args) ",") ") "
(expr->Octave body #t) ")"))
(V:Op (name args)
(let* ((infix? (case name
((+ - * / >= > < <=) #t)
(else #f)))
(op name))
(cond ((null? args)
(case name
((NONE) (list name))
(else (list name "()"))))
((null? (cdr args))
(list op "(" (value->Octave (car args)) ")"))
((and infix? (null? (cddr args)))
(list "(" (value->Octave (car args)) ")" op "(" (value->Octave (cadr args)) ")"))
(infix?
(let ((op (case op ((+) 'sum) ((*) 'prod) (else op))))
(list op "([" (intersperse (map value->Octave args) ",") "])")))
(else
(list op "(" (intersperse (map value->Octave args) ",") ")")))))
(V:Ifv (test ift iff)
(list "(ifv (" (value->Octave test) ", "
(value->Octave ift ) ", "
(value->Octave iff) "))"))
))
(define (prelude/Octave #!key (solver 'lsode))
`(
#<