source: project/release/4/9ML-toolkit/trunk/ivp-octave.scm @ 25514

Last change on this file since 25514 was 25514, checked in by Ivan Raikov, 9 years ago

9ML-toolkit: support for lsode solver in octave

File size: 4.3 KB
Line 
1;;
2;;  NineML IVP code generator for Octave.
3;;
4;;
5;; Copyright 2010-2011 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-octave
24
25        (ivp-octave)
26
27        (import scheme chicken )
28
29(import (only extras sprintf)
30        (only files make-pathname pathname-directory pathname-file absolute-pathname? )
31        (only data-structures conc alist-ref intersperse ->string)
32        (only posix current-directory)
33        (only srfi-1 filter list-index)
34        (only srfi-13 string-concatenate))
35(require-extension setup-api datatype signal-diagram 9ML-repr )
36
37
38(define nl "\n")
39
40
41(define (octave-m ivp-id start ivar hvar dvars events ic imax log-path adir solver-path)
42
43  (let* ((states (cons ivar dvars))
44         (parameters (filter (lambda (x) (not (member (car x) states))) ic)))
45
46  `(,(sprintf "autoload (~S,~S);~%" (->string ivp-id) (make-pathname adir solver-path))
47    ,(sprintf "~A;~%~%" ivp-id)
48   
49    ,(sprintf "~A_initial = struct(~A);~%~%" ivp-id
50              (string-concatenate
51               (intersperse
52                (map (lambda (x)
53                       (let ((n (car x)) (v (cdr x)))
54                         (let ((v (if (boolean? v) (if v "true" "false") v)))
55                           (sprintf "\"~A\",~A" n v))))
56                     (cons (cons ivar start) ic))
57                ",")))
58
59    ,(sprintf "i = 1; t = 0; tmax = ~A;~%" imax)
60    ,(sprintf "~A_state = ~A_initial;~%" ivp-id ivp-id)
61    ,(sprintf "data = [[~A]];~%" 
62              (string-concatenate (intersperse (map (lambda (s) (sprintf "~A_state.~A" ivp-id s)) states) ",") ))
63
64    ,(sprintf "tic;~%")
65    ,(sprintf "while (t < tmax) ~%")
66
67    ,(sprintf "  y = ~A(~A_state);~%" ivp-id ivp-id)
68
69    ,(if (null? events) '()
70         `(,(sprintf "  if (~A)~%" (intersperse (map (lambda (x) (sprintf "y.~A" x)) events) " || "))
71           ,(sprintf "    ~A_state.~A = 1e-4;~%" ivp-id hvar)
72           ,(sprintf "    y = ~A(~A_state);~%" ivp-id ivp-id)
73           ,(sprintf "  else~%")
74           ,(sprintf "    if (~A_state.~A < 0.25)~%" ivp-id hvar)
75           ,(sprintf "      ~A_state.~A = ~A_state.~A+1e-2;~%" ivp-id hvar ivp-id hvar)
76           ,(sprintf "    endif~%")
77           ,(sprintf "  endif~%")
78           ))
79
80    ,(map (lambda (s) (sprintf "  ~A_state.~A = y.~A;~%" ivp-id s s)) states) 
81    ,(sprintf "  data = [data; [~A]];~%" 
82              (string-concatenate (intersperse (map (lambda (s) (sprintf "y.~A" s)) states) ",") ))
83   
84    ,(sprintf "  t  = y.~A;~%" ivar )
85    ,(sprintf "  i = i+1;~%")
86    ,(sprintf "endwhile~%")
87    ,(sprintf "toc;~%")
88
89    ,(sprintf "save (\"-ascii\", ~S, \"data\");~%" log-path)
90    )
91  ))
92
93
94
95
96(define (ivp-octave prefix ivp-id hvar ivar dvars pvars events start end ic sd)
97
98  (let* (                       
99         (N (+ 1 (length dvars)))
100         (P (- (length ic) N))
101         (dir        (or (pathname-directory prefix) "."))
102         (adir       (if (absolute-pathname? dir) dir
103                         (make-pathname (current-directory) dir)))
104         (shared-dir         (chicken-home))
105         (solver-path        (make-pathname dir (sprintf "~A_solver.m" ivp-id)))
106         (m-path             (make-pathname dir (sprintf "~A_~A.m" (pathname-file prefix) ivp-id)))
107         (log-path           (make-pathname dir (sprintf "~A_~A.log" (pathname-file prefix) ivp-id)))
108         (octave-path        "octave")
109         )
110   
111    (make (
112           (solver-path (prefix)
113                        (with-output-to-file solver-path
114                          (lambda () 
115                            (codegen/Octave ivp-id sd solver: 'lsode))) )
116           
117           (m-path ()
118                   (with-output-to-file m-path
119                     (lambda () 
120                       (print-fragments
121                        `(,(octave-m ivp-id start ivar hvar dvars events ic end log-path adir solver-path))))) )
122           
123           
124           (log-path (m-path solver-path)
125                     (run (octave -q ,m-path)) )
126
127#|         
128           (clean ()
129                  (run rm ,solver-path ,clib-path ,h-path
130                       ,mlb-path ,clib-so-path ,lib-so-path
131                       ,oct-cc-path ,oct-run-path ,oct-initial-path
132                       ,oct-open-path ,oct-close-path ,m-path ,log-path))
133|#
134           
135           )
136      (list solver-path m-path log-path) )))
137
138)
139
Note: See TracBrowser for help on using the repository browser.