Commit f844645d authored by liam's avatar liam
Browse files

More definitions for ordinary differential equations, and example.

Allocating stepper causes memory fault.


git-svn-id: svn+ssh://pop/opt/space/mathematics/gsl/trunk@3234 a3d8a0fb-c1db-0310-ace7-a616afeb9e30
parent 46d28c60
;********************************************************
; file: evolution.lisp
; description: Evolution functions for ODE integration.
; date: Sun Sep 30 2007 - 14:31
; author: Liam Healy
; modified: Sun Sep 30 2007 - 15:50
;********************************************************
;;; $Id: $
(in-package :gsl)
(defun-gsl allocate-evolution (dimension)
"gsl_odeiv_evolve_alloc"
((dimension :size))
:c-return :pointer
:documentation
"Allocate a new instance of an evolution function
for a system of dimension dimensions and return the pointer.")
(defun-gsl apply-evolution
(evolve control step dydt time max-time step-size y)
"gsl_odeiv_evolve_apply"
((evolve :pointer) (control :pointer) (step :pointer)
(dydt :pointer) (time :double) (max-time :double)
(step-size :double) (y :pointer))
:documentation
"Advance the system (@var{e}, @var{dydt}) from time
@var{t} and position @var{y} using the stepping function @var{step}.
The new time and position are stored in @var{time} and @var{y} on output.
The initial step-size is taken as @var{step-size}, but this will be modified
using the control function @var{c} to achieve the appropriate error
bound if necessary. The routine may make several calls to @var{step} in
order to determine the optimum step-size. If the step-size has been
changed the value of @var{h} will be modified on output. The maximum
time @var{max-time} is guaranteed not to be exceeded by the time-step. On the
final time-step the value of @var{time} will be set to @var{t1} exactly.")
(defun-gsl reset-evolution (evolve)
"gsl_odeiv_evolve_reset"
((evolve :pointer))
:documentation
"Reset the evolution function @var{evolve}. It should be used
whenever the next use of @var{evolve} will not be a continuation of a
previous step.")
(defun-gsl free-evolution (evolve)
"gsl_odeiv_evolve_free"
((evolve :pointer))
:c-return :void
:documentation
"Frees all the memory associated with the evolution function.")
......@@ -3,7 +3,7 @@
; description: Example ODE
; date: Sat Sep 29 2007 - 17:49
; author: Liam Healy
; modified: Sat Sep 29 2007 - 18:51
; modified: Sun Sep 30 2007 - 16:58
;********************************************************
;;; $Id: $
......@@ -28,6 +28,34 @@
(double-to-cl dfdy 3)
(* -1 mu (- (* (double-to-cl y 0) (double-to-cl y 0)) 1.0d0))))
(let ((stepper (step-allocate *step-rk8pd* 2)
(step-free stepper))))
(defun integrate-vanderpol (max-time &optional (step-size 1.0d-6))
(let ((stepper (step-allocate *step-rk8pd* 2))
(control (new-y-control 1.0d-6 0.0d0))
(evolve (allocate-evolution 2))
(mu 10.0d0))
(declare (special mu))
(cffi:with-foreign-objects
((vanderpol 'ode-system) (dependent :double 2))
(setf
(cffi:foreign-slot-value vanderpol 'ode-system 'function)
(cffi:callback vanderpol-function)
(cffi:foreign-slot-value vanderpol 'ode-system 'jacobian)
(cffi:callback vanderpol-jacobian)
(cffi:foreign-slot-value vanderpol 'ode-system 'dimension)
2
(cffi:foreign-slot-value vanderpol 'ode-system 'parameters)
(cffi:null-pointer)
(double-to-cl dependent 0) 1.0d0
(double-to-cl dependent 1) 0.0d0)
(loop for time from 0.0d0 below max-time
do
(apply-evolution
evolve control stepper vanderpol
time max-time step-size dependent)
(format t "~&~d ~d ~d"
time
(double-to-cl dependent 0)
(double-to-cl dependent 1))))
(free-evolution evolve)
(free-control control)
(step-free stepper)))
......@@ -3,7 +3,7 @@
; description: ODE system setup
; date: Sun Apr 15 2007 - 14:19
; author: Liam Healy
; modified: Sat Sep 29 2007 - 17:54
; modified: Sun Sep 30 2007 - 15:47
;********************************************************
;;; $Id: $
......@@ -17,7 +17,7 @@
(dimension :size)
(parameters :pointer))
(export 'def-ode-function 'def-jacobian-function)
(export '(def-ode-function def-jacobian-function))
(defmacro def-ode-function (name (time dependent derivatives) &body body)
"Define a function that will evaluate the right-hand sides (derivatives)
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment