Wiki
Download
Manual
Eggs
API
Tests
Bugs
show
edit
history
You can edit this page using
wiki syntax
for markup.
Article contents:
== Outdated egg! This is an egg for CHICKEN 4, the unsupported old release. You're almost certainly looking for [[/eggref/5/seulex|the CHICKEN 5 version of this egg]], if it exists. If it does not exist, there may be equivalent functionality provided by another egg; have a look at the [[https://wiki.call-cc.org/chicken-projects/egg-index-5.html|egg index]]. Otherwise, please consider porting this egg to the current version of CHICKEN. [[tags: eggs]] [[toc:]] == seulex === Description The Chicken {{seulex}} library provides bindings to the Fortran {{SEULEX}} solver for systems of stiff differential and differential algebraic equations of the form {{MY'=F(X,Y)}}. SEULEX was written by E. Hairer AND G. Wanner and stands for Semi-implicit EULer EXtrapolation and computes its error estimates by means of polynomial extrapolation. Extrapolation solvers can achieve good performance for certain types of problems where high precision is required. === Installation notes The Chicken {{seulex}} library must be compiled along with the {{SEULEX}} Fortran source code. Therefore, the user must have a Fortran compiler such as gfortran or g95 installed on their system. The following environment variables can be used to control the Fortran compilation process: ; {{FORTRAN_COMPILER}} : path to Fortran compiler (default is "gfortran") ; {{FORTRAN_FLAGS}} : flags to be passed to the Fortran compiler (default is "-fno-automatic -fPIC -g") ; {{FORTRAN_LIBS}} : Fortran libraries to link to (default is "-lgfortran -lblas -llapack") === Library procedures <procedure>(seulex-create-solver XINIT YINIT FCN [JACOBIAN] [AUTONOMOUS] [USER-DATA] [DENSITY-COMPONENTS] [RELTOL] [ABSTOL]) => SEULEX-SOLVER</procedure> Creates and initializes an object representing a problem to be solved with the SEULEX solver. Arguments {{XINIT}} and {{YINIT}} represent the initial conditions: {{XINIT}} is a scalar value for the independent variable, and {{YINIT}} is an SRFI-4 {{f64vector}} value containing the initial dependent variable values. Argument {{FCN}} is used to compute the right-hand side function {{F}} and must be a procedure of the following form: (LAMBDA T YY DATA) or (LAMBDA T YY) depending on whether the {{USER-DATA}} optional argument is set, where ; {{T}} : real-valued independent variable ; {{YY}} : SRFI-4 {{f64vector}} with current variable values ; {{DATA}} : is a user data object (if set) This procedure must return a SRFI-4 {{f64vector}} containing the derivative values. Optional keyword argument {{JACOBIAN}} is a procedure which will be used to compute the partial derivatives of {{F(X,Y)}} with respect to {{Y}}. If not given, it is computed internally by finite differences. Optional keyword argument {{AUTONOMOUS}} is a boolean value that indicates whether {{F(X,Y)}} is independent of {{X}} (autonomous) or not (non-autonomous). The default is {{#t}} (autonomous). Optional keyword argument {{USER-DATA}} is an object that will be passed as an additional argument to the right-hand side function. Optional keyword argument {{DENSITY-COMPONENTS}} must be an SRFI-4 {{s32vector}} which indicates the components for which dense output is required. Optional keyword arguments {{RELTOL}} and {{ABSTOL}} specify relative and absolute error tolerance, respectively. These both default to 1e-4. <procedure>(seulex-destroy-solver SEULEX-SOLVER)</procedure> Deallocates the memory associated with the given solver. <procedure>(seulex-solve SEULEX-SOLVER T)</procedure> Integrates the system over an interval in the independent variable. This procedure returns either when the given {{T}} is reached, or when a root is found. <procedure>(seulex-yy SEULEX-SOLVER)</procedure> Returns the SRFI-4 {{f64vector}} value of current state values of the system. <procedure>(seulex-nfcn SEULEX-SOLVER)</procedure> Returns the number of function evaluations done so far. <procedure>(seulex-njac SEULEX-SOLVER)</procedure> Returns the number of Jacobian function evaluations done so far. <procedure>(seulex-nstep SEULEX-SOLVER)</procedure> Returns the number of computed steps. <procedure>(seulex-ndec SEULEX-SOLVER)</procedure> Returns the number of LU decompositions. <procedure>(seulex-nsol SEULEX-SOLVER)</procedure> Returns the number of backward-forward substitutions. === Example ;; ;; Hodgkin-Huxley model ;; (use mathh seulex srfi-4) (define neg -) (define pow expt) (define TEND 500.0) ;; Model parameters (define (I_stim t) 10) (define C_m 1) (define E_Na 50) (define E_K -77) (define E_L -54.4) (define gbar_Na 120) (define gbar_K 36) (define g_L 0.3) ;; Rate functions (define (amf v) (* 0.1 (/ (+ v 40) (- 1.0 (exp (/ (neg (+ v 40)) 10)))))) (define (bmf v) (* 4.0 (exp (/ (neg (+ v 65)) 18)))) (define (ahf v) (* 0.07 (exp (/ (neg (+ v 65)) 20)))) (define (bhf v) (/ 1.0 (+ 1.0 (exp (/ (neg (+ v 35)) 10))))) (define (anf v) (* 0.01 (/ (+ v 55) (- 1 (exp (/ (neg (+ v 55)) 10)))))) (define (bnf v) (* 0.125 (exp (/ (neg (+ v 65)) 80)))) ;; State functions (define (minf v) (* 0.5 (+ 1 (tanh (/ (- v v1) v2))))) (define (winf v) (* 0.5 (+ 1 (tanh (/ (- v v3) v4))))) (define (lamw v) (* phi (cosh (/ (- v v3) (* 2 v4))))) ;; Model equations (define (rhs t yy) (let ((v (f64vector-ref yy 0)) (m (f64vector-ref yy 1)) (h (f64vector-ref yy 2)) (n (f64vector-ref yy 3))) ;; transition rates at current step (let ((am (amf v)) (an (anf v)) (ah (ahf v)) (bm (bmf v)) (bn (bnf v)) (bh (bhf v)) (g_Na (* gbar_Na (* h (pow m 3)))) (g_K (* gbar_K (pow n 4)))) (let ( ;; currents (I_Na (* (- v E_Na) g_Na)) (I_K (* (- v E_K) g_K)) (I_L (* g_L (- v E_L)))) (let ( ;; state equations (dm (- (* am (- 1 m)) (* bm m))) (dh (- (* ah (- 1 h)) (* bh h))) (dn (- (* an (- 1 n)) (* bn n))) (dv (/ (- (I_stim t) I_L I_Na I_K) C_m)) ) (f64vector dv dm dh dn) ))) )) (let ((yy (f64vector -65 0.052 0.596 0.317)) ;; v m h n ;; Integration limits (t0 0.0) (tf TEND) (dt 1e-2)) ;; solver initialization (let ((solver (seulex-create-solver t0 yy rhs abstol: 1e-4 reltol: 1e-4))) ;; In loop, call solver, print results, and test for error. (let recur ((tnext (+ t0 dt)) (iout 1)) (let ((flag (seulex-solve solver tnext))) (if (negative? flag) (error 'main "SEULEX solver error" flag)) (print-results solver tnext) (if (< tnext tf) (recur (+ tnext dt) (+ 1 iout))) )) (seulex-destroy-solver solver) (define (print-results solver t) (let ((yy (seulex-yy solver))) (printf "~A ~A ~A ~A ~A~%" t (f64vector-ref yy 0) (f64vector-ref yy 1) (f64vector-ref yy 2) (f64vector-ref yy 3) ))) === Version History * 1.0 Initial release === License The SEULEX solver was written by E. HAIRER AND G. WANNER, UNIVERSITE DE GENEVE, DEPT. DE MATHEMATIQUES Chicken Scheme bindings for SEULEX are copyright 2011-2012 Ivan Raikov. 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 <http://www.gnu.org/licenses/>.
Description of your changes:
I would like to authenticate
Authentication
Username:
Password:
Spam control
What do you get when you multiply 7 by 2?