Commit 3f6e03fa authored by liam's avatar liam
Browse files

ODE integration now done, example works. Pointers rather than the

doubles themselves passed to #'apply-evolution.  Interface could use
work, at least a helper macro.


git-svn-id: svn+ssh://pop/opt/space/mathematics/gsl/trunk@3238 a3d8a0fb-c1db-0310-ace7-a616afeb9e30
parent 0275853d
...@@ -3,7 +3,7 @@ ...@@ -3,7 +3,7 @@
; description: Definition of GSLL system ; description: Definition of GSLL system
; date: ; date:
; author: Liam Healy ; author: Liam Healy
; modified: Sun Feb 11 2007 - 16:14 ; modified: Thu Oct 18 2007 - 22:18
;******************************************************** ;********************************************************
;;; $Id: $ ;;; $Id: $
...@@ -151,4 +151,12 @@ ...@@ -151,4 +151,12 @@
(:file "read-write" :depends-on (histogram)) (:file "read-write" :depends-on (histogram))
(:file "probability-distribution" :depends-on (histogram)) (:file "probability-distribution" :depends-on (histogram))
(:file "ntuple"))) (:file "ntuple")))
(:file "monte-carlo" :depends-on (init data random numerical-integration)))) (:file "monte-carlo" :depends-on (init data random numerical-integration))
(:module ordinary-differential-equations
:depends-on (init)
:components
((:file "ode-system")
(:file "stepping")
(:file "control")
(:file "evolution")
(:file "ode-example" :depends-on (ode-system stepping))))))
...@@ -268,7 +268,7 @@ svn checkout svn://common-lisp.net/project/gsll/subversion/trunk</pre> ...@@ -268,7 +268,7 @@ svn checkout svn://common-lisp.net/project/gsll/subversion/trunk</pre>
<tr> <tr>
<td><a href="Ordinary-Differential-Equations.html">Ordinary <td><a href="Ordinary-Differential-Equations.html">Ordinary
Differential Equations</a></td> Differential Equations</a></td>
<td>Not done</td> <td>Done</td>
</tr> </tr>
<tr> <tr>
<td><a href="Interpolation.html">Interpolation</a></td> <td><a href="Interpolation.html">Interpolation</a></td>
......
...@@ -3,7 +3,7 @@ ...@@ -3,7 +3,7 @@
; description: Evolution functions for ODE integration. ; description: Evolution functions for ODE integration.
; date: Sun Sep 30 2007 - 14:31 ; date: Sun Sep 30 2007 - 14:31
; author: Liam Healy ; author: Liam Healy
; modified: Sun Sep 30 2007 - 15:50 ; modified: Thu Oct 18 2007 - 22:18
;******************************************************** ;********************************************************
;;; $Id: $ ;;; $Id: $
...@@ -21,8 +21,8 @@ ...@@ -21,8 +21,8 @@
(evolve control step dydt time max-time step-size y) (evolve control step dydt time max-time step-size y)
"gsl_odeiv_evolve_apply" "gsl_odeiv_evolve_apply"
((evolve :pointer) (control :pointer) (step :pointer) ((evolve :pointer) (control :pointer) (step :pointer)
(dydt :pointer) (time :double) (max-time :double) (dydt :pointer) (time :pointer) (max-time :double)
(step-size :double) (y :pointer)) (step-size :pointer) (y :pointer))
:documentation :documentation
"Advance the system (@var{e}, @var{dydt}) from time "Advance the system (@var{e}, @var{dydt}) from time
@var{t} and position @var{y} using the stepping function @var{step}. @var{t} and position @var{y} using the stepping function @var{step}.
......
...@@ -3,10 +3,14 @@ ...@@ -3,10 +3,14 @@
; description: Example ODE ; description: Example ODE
; date: Sat Sep 29 2007 - 17:49 ; date: Sat Sep 29 2007 - 17:49
; author: Liam Healy ; author: Liam Healy
; modified: Sun Sep 30 2007 - 16:58 ; modified: Thu Oct 18 2007 - 22:39
;******************************************************** ;********************************************************
;;; $Id: $ ;;; $Id: $
;;; Interfaces to GSL are still somewhat crude, as can be seen in the
;;; van der Pol example implemented below.
;;; Thu Oct 18 2007
(in-package :gsl) (in-package :gsl)
(def-ode-function vanderpol-function (time y dydt) (def-ode-function vanderpol-function (time y dydt)
...@@ -28,14 +32,21 @@ ...@@ -28,14 +32,21 @@
(double-to-cl dfdy 3) (double-to-cl dfdy 3)
(* -1 mu (- (* (double-to-cl y 0) (double-to-cl y 0)) 1.0d0)))) (* -1 mu (- (* (double-to-cl y 0) (double-to-cl y 0)) 1.0d0))))
(defparameter *max-iter* 2000)
(defun integrate-vanderpol (max-time &optional (step-size 1.0d-6)) (defun integrate-vanderpol (max-time &optional (step-size 1.0d-6))
"Integrate the van der Pol oscillator as given in Section 25.5 of the
GSL manual. To reproduce that example, (integrate-vanderpol 100.0d0)."
(let ((stepper (step-allocate *step-rk8pd* 2)) (let ((stepper (step-allocate *step-rk8pd* 2))
(control (new-y-control 1.0d-6 0.0d0)) (control (new-y-control 1.0d-6 0.0d0))
(evolve (allocate-evolution 2)) (evolve (allocate-evolution 2))
(mu 10.0d0)) (mu 10.0d0)
(time 0.0d0)
(iter 0))
(declare (special mu)) (declare (special mu))
(cffi:with-foreign-objects (cffi:with-foreign-objects
((vanderpol 'ode-system) (dependent :double 2)) ((vanderpol 'ode-system) (dependent :double 2)
(ctime :double) (cstep :double))
(setf (setf
(cffi:foreign-slot-value vanderpol 'ode-system 'function) (cffi:foreign-slot-value vanderpol 'ode-system 'function)
(cffi:callback vanderpol-function) (cffi:callback vanderpol-function)
...@@ -46,13 +57,16 @@ ...@@ -46,13 +57,16 @@
(cffi:foreign-slot-value vanderpol 'ode-system 'parameters) (cffi:foreign-slot-value vanderpol 'ode-system 'parameters)
(cffi:null-pointer) (cffi:null-pointer)
(double-to-cl dependent 0) 1.0d0 (double-to-cl dependent 0) 1.0d0
(double-to-cl dependent 1) 0.0d0) (double-to-cl dependent 1) 0.0d0
(loop for time from 0.0d0 below max-time (double-to-cl ctime) time
do (double-to-cl cstep) step-size)
(loop (when (or (>= time max-time) (> iter *max-iter*)) (return))
(apply-evolution (apply-evolution
evolve control stepper vanderpol evolve control stepper vanderpol
time max-time step-size dependent) ctime max-time cstep dependent)
(format t "~&~d ~d ~d" (setf time (double-to-cl ctime))
(incf iter)
(format t "~&~12,6f~10t~12,6f~24t~12,6f"
time time
(double-to-cl dependent 0) (double-to-cl dependent 0)
(double-to-cl dependent 1)))) (double-to-cl dependent 1))))
......
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