source: project/release/4/9ML-toolkit/trunk/ivp-lib.scm @ 29952

Last change on this file since 29952 was 29952, checked in by Ivan Raikov, 8 years ago

9ML-toolkit: bringing octave and scheme backends up to date; changed Izhikevich FS example to use heaviside function

File size: 9.6 KB
Line 
1;;
2;;  An IVP solver for NineML.
3;;
4;;
5;; Copyright 2010-2013 Ivan Raikov and the Okinawa Institute of
6;; Science and Technology.
7;;
8;; This program is free software: you can redistribute it and/or
9;; modify it under the terms of the GNU General Public License as
10;; published by the Free Software Foundation, either version 3 of the
11;; License, or (at your option) any later version.
12;;
13;; This program is distributed in the hope that it will be useful, but
14;; WITHOUT ANY WARRANTY; without even the implied warranty of
15;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
16;; General Public License for more details.
17;;
18;; A full copy of the GPL license can be found at
19;; <http://www.gnu.org/licenses/>.
20;;
21
22
23(module 9ML-ivp-lib
24
25        (
26         ivp-verbose  simulation-platform 
27         make-ivp-cgen-hook make-ivp-data-hook make-ivp-plot-hook
28         )
29
30(import scheme chicken)
31
32(require-extension datatype matchable static-modules miniML miniMLsyntax miniMLvalue miniMLeval)
33(require-extension object-graph signal-diagram)
34(require-extension 9ML-parse 9ML-eval 9ML-plot)
35(require-extension 9ML-ivp-octave 9ML-ivp-chicken 9ML-ivp-mlton 9ML-ivp-octave-mlton )
36
37(import (only files pathname-file)
38        (only utils system*)
39        (only extras fprintf)
40        (only posix process-fork)
41        (only data-structures string-intersperse ->string alist-ref)
42        (only srfi-1 lset-intersection lset-difference)
43        )
44
45(define lookup-def 
46  (lambda (k lst . rest)
47    (let-optionals rest ((default #f))
48      (alist-ref k lst eq? default))))
49
50
51(define simulation-platform (make-parameter #f))
52(define ivp-verbose (make-parameter 0))
53
54(define (d fstr . args)
55  (let ([port (current-error-port)])
56    (if (positive? (ivp-verbose)) 
57        (begin (apply fprintf port fstr args)
58               (flush-output port) ) )))
59
60
61(define (run:execute explist)
62  (define (smooth lst)
63    (let ((slst (map ->string lst)))
64      (string-intersperse (cons (car slst) (cdr slst)) " ")))
65  (for-each (lambda (cmd) (system (->string cmd)))
66            (map smooth explist)))
67
68
69(define (run:execute* explist)
70  (define (smooth lst)
71    (let ((slst (map ->string lst)))
72      (string-intersperse (cons (car slst) (cdr slst)) " ")))
73  (for-each (lambda (cmd) (system* "~a" cmd))
74            (map smooth explist)))
75
76
77(define-syntax run
78  (syntax-rules ()
79    ((_ exp ...)
80     (begin
81       (d "running ~A ...~%" (list `exp ...))
82       (run:execute* (list `exp ...))))))
83
84
85(define-syntax run-
86  (syntax-rules ()
87    ((_ exp ...)
88     (begin
89       (d "running ~A ...~%" (list `exp ...))
90       (run:execute (list `exp ...))))))
91
92
93;; Use args:usage to generate a formatted list of options (from OPTS),
94;; suitable for embedding into help text.
95
96(define nl "\n")
97
98
99(define (construct-ivp prefix name sxml-tuple)
100
101  (let ((sexpr 
102         (let ((sexpr (sxml-value->sexpr sxml-tuple)))
103           (case (car sexpr)
104             ((ivp)
105              (and (pair? (cdr sexpr))
106                   (case (cadr sexpr)
107                     ((initial)  (cddr sexpr))
108                     (else #f))))
109             (else #f)))))
110
111      (and sexpr
112           (match-let (((ivar hvar) (cdr sexpr)))
113
114             (let* (
115                    (diagram+initial (sexpr->diagram+initial hvar (car sexpr)))
116                    (sd (construct (car diagram+initial)))
117                    (ic (cadr diagram+initial))
118                    (fields (caddr diagram+initial))
119                    )
120
121               (if (not (alist-ref ivar ic))
122                   (error 'construct-ivp "IVP independent variable is not present in given system" name ivar))
123               (if (not (alist-ref hvar ic))
124                   (error 'construct-ivp "IVP step variable is not present in given system" name hvar))
125
126               (let* ((dfe (dataflow sd '()))
127                      (dvars (lset-difference eq?
128                                              (lset-intersection eq? (alist-ref 'in dfe) (alist-ref 'out dfe))
129                                              (list ivar)))
130                      (pvars (lset-difference eq? (alist-ref 'in dfe) (cons ivar dvars)))
131                      (events (reverse (events sd)))
132                      )
133
134                 `(
135                   (signal-diagram . ,sd)
136                   (initial-conditions . ,ic)
137                   (dfe    . ,dfe)
138                   (ivar   . ,ivar)
139                   (hvar   . ,hvar)
140                   (dvars  . ,dvars)
141                   (pvars  . ,pvars)
142                   (events . ,events)
143                   (fields . ,fields)
144                   )
145                 ))
146             ))
147      ))
148                   
149                 
150
151(define (generate-ivp-code prefix ivp-id sxml-tuple #!key (platform 'chicken) (method 'rk3))
152
153  (let ((sdinfo (construct-ivp prefix ivp-id sxml-tuple)))
154
155    (let (
156          (sd      (lookup-def 'signal-diagram sdinfo))
157          (ic      (lookup-def 'initial-conditions sdinfo))
158          (ivar    (lookup-def 'ivar sdinfo))
159          (hvar    (lookup-def 'hvar sdinfo))
160          (dvars   (lookup-def 'dvars sdinfo))
161          (pvars   (lookup-def 'pvars sdinfo))
162          (events  (lookup-def 'events sdinfo))
163          )
164
165      (d "generate-ivp-code: ic = ~A~%" ic)
166     
167      (case platform
168       
169        ((octave)
170         (begin
171           (ivp-octave-codegen prefix ivp-id hvar ivar dvars pvars events ic sd)
172           `((ivp-id . ,ivp-id) . ,sdinfo)))
173       
174        ((octave/mlton octave-mlton)
175         (begin
176           (ivp-octave-mlton-codegen prefix ivp-id hvar ivar dvars pvars ic sd)
177           `((ivp-id . ,ivp-id) . ,sdinfo)))
178       
179        ((mlton)
180         (begin
181           (ivp-mlton-codegen prefix ivp-id ivar dvars pvars ic sd method)
182           `((ivp-id . ,ivp-id) . ,sdinfo)))
183       
184        ((chicken)
185         (begin
186           (ivp-chicken-codegen prefix ivp-id ivar dvars pvars events ic sd method)
187           `((ivp-id . ,ivp-id) . ,sdinfo)))
188       
189        (else (error 'generate-ivp-code "unknown platform"  platform))
190       
191        ))
192    ))
193   
194                 
195
196(define (generate-ivp-table prefix ivp-id sxml-tuple #!key (platform 'chicken))
197
198               
199  (let ((sexpr 
200         (let ((sexpr (sxml-value->sexpr sxml-tuple)))
201           (case (car sexpr)
202             ((ivp)
203              (and (pair? (cdr sexpr))
204                   (case (cadr sexpr)
205                     ((run)  (cddr sexpr))
206                     (else #f))))
207             (else #f)))))
208
209    (d "generate-ivp-table: sexpr = ~A~%" (sxml-value->sexpr sxml-tuple))
210
211      (and sexpr
212           (match-let (((ivar hvar start end) (cdr sexpr)))
213
214             (let* ((diagram+initial (sexpr->diagram+initial hvar (car sexpr)))
215                    (sd (construct (car diagram+initial)))
216                    (ic (cadr diagram+initial)))
217
218               (d "generate-ivp-table: ic = ~A~%" ic)
219
220               (if (not (alist-ref ivar ic))
221                   (error 'generate-ivp-table "IVP independent variable is not present in given system" ivar))
222               (if (not (alist-ref hvar ic))
223                   (error 'generate-ivp-table "IVP step variable is not present in given system" hvar))
224               (let* ((dfe (dataflow sd '()))
225                      (dvars (lset-difference eq?
226                                              (lset-intersection eq? (alist-ref 'in dfe) (alist-ref 'out dfe))
227                                              (list ivar)))
228                      (pvars (lset-difference eq? (alist-ref 'in dfe) (cons ivar dvars)))
229                      (events (events sd))
230                      )
231
232
233                 (case platform
234                   
235                   ((octave)
236                    (process-fork
237                     (lambda () (ivp-octave prefix ivp-id hvar ivar dvars pvars events start end ic sd)))
238                    (list ivp-id ivar dvars) )
239
240                   ((octave/mlton octave-mlton)
241                    (process-fork
242                     (lambda () (ivp-octave-mlton prefix ivp-id hvar ivar dvars pvars start end ic sd)))
243                    (list ivp-id ivar dvars) )
244                   
245                   ((mlton)
246                    (process-fork
247                     (lambda () (ivp-mlton  prefix ivp-id ivar dvars pvars start end ic sd 'rk3)))
248                    (list ivp-id ivar dvars) )
249
250                   ((chicken)
251                    (process-fork
252                     (lambda () (ivp-chicken prefix ivp-id ivar dvars pvars events start end ic sd)))
253                    (list ivp-id ivar dvars) )
254
255                   (else (error 'generate-ivp-table "unknown platform"  platform))
256                   
257                   )))
258             )))
259  )
260
261
262
263(define (make-ivp-cgen-hook ivp-env)
264  (lambda (prefix name label value)
265
266    (cond
267
268     ((or (and (string? label) (string=? label "ivp"))
269          (and (pair? label) (string=? (car label) "ivp"))) ;; value is an IVP
270
271      (let ((ivp-info (generate-ivp-code prefix name value method: 'rkfe platform: (simulation-platform))))
272        (ivp-env (cons `(,(string->symbol name) . ,ivp-info) (ivp-env)))
273        ivp-info
274        ))
275     
276     (else #f)
277     ))
278  )
279           
280
281(define (make-ivp-data-hook #!key (ivp #f) (diagram #f))
282  (lambda (prefix name label value)
283
284    (cond
285     ((and diagram
286           (or (and (string? label) (string=? label "diagram")) 
287               (and (pair? label) (string=? (car label) "diagram")))) ;; value is a diagram
288      (let* ((diagram-id (gensym 'diagram))
289             (diagram-link `(img (@ (src ,(sprintf "~A.png" diagram-id))) (alt "NineML diagram"))))
290        (plot-diagram prefix diagram-id value)
291        diagram-link
292        ))
293     
294     ((and ivp (or (and (string? label) (string=? label "ivp"))
295                   (and (pair? label) (string=? (car label) "ivp")))) ;; value is an IVP
296      (let ((ivp-id (gensym (string->symbol (string-append (->string name) "ivp")))))
297        (let ((ivp-info (generate-ivp-table prefix ivp-id value platform: (simulation-platform))))
298          ivp-info
299        ))
300      )
301     
302     (else #f))))
303           
304
305(define (make-ivp-plot-hook #!key (ivp #f) (plot-format 'png) (plot-index #f))
306  (lambda (prefix name label value)
307    (cond
308
309     ((and ivp (or (and (string? label) (string=? label "ivp"))
310                   (and (pair? label) (string=? (car label) "ivp")))) ;; value is an IVP
311      (let* ((ivp-id (gensym 'ivp))
312             (ivp-plot-link `(img (@ (src ,(sprintf "~A_~A.png" (pathname-file prefix) ivp-id)) (alt "NineML IVP plot")))))
313
314        (let ((ivp-info (generate-ivp-table prefix ivp-id value platform: (simulation-platform))))
315          (plot-ivp prefix ivp-info format: plot-format index: plot-index)
316          ivp-plot-link)
317        ))
318     
319     (else #f))))
320           
321
322
323)
Note: See TracBrowser for help on using the repository browser.