diff --git a/ordinary-differential-equations/evolution.lisp b/ordinary-differential-equations/evolution.lisp
index aa9e69725ff79c49bdbf528de71da41c87e7a585..49839b6ac422c0a08e1a85c44cb8c3dc3605b19e 100644
--- a/ordinary-differential-equations/evolution.lisp
+++ b/ordinary-differential-equations/evolution.lisp
@@ -1,6 +1,6 @@
 ;; Evolution functions for ODE integration.
 ;; Liam Healy, Sun Sep 30 2007 - 14:31
-;; Time-stamp: <2009-02-15 08:35:23EST evolution.lisp>
+;; Time-stamp: <2009-02-15 08:58:47EST evolution.lisp>
 ;; $Id$
 
 (in-package :gsl)
@@ -13,11 +13,11 @@
   :initialize-args nil)
 
 (defmfun apply-evolution
-    (evolve control step time max-time step-size y)
+    (evolution time y step-size control stepper max-time)
   "gsl_odeiv_evolve_apply"
-  (((mpointer evolve) :pointer) ((mpointer control) :pointer)
-   ((mpointer step) :pointer)
-   ((callback-struct step) :pointer) ((c-pointer time) :pointer)
+  (((mpointer evolution) :pointer) ((mpointer control) :pointer)
+   ((mpointer stepper) :pointer)
+   ((callback-struct stepper) :pointer) ((c-pointer time) :pointer)
    (max-time :double)
    ((c-pointer step-size) :pointer) ((c-pointer y) :pointer))
   :inputs (time step-size y)
diff --git a/ordinary-differential-equations/ode-example.lisp b/ordinary-differential-equations/ode-example.lisp
index a00f258445c48dcf405d4eac26ad370a22a31103..5aa1819072aa4cc75555a89d69e3a4d3c79921ee 100644
--- a/ordinary-differential-equations/ode-example.lisp
+++ b/ordinary-differential-equations/ode-example.lisp
@@ -1,6 +1,6 @@
 ;; Example ODE                               
 ;; Liam Healy Sat Sep 29 2007 - 17:49
-;; Time-stamp: <2009-02-15 08:39:44EST ode-example.lisp>
+;; Time-stamp: <2009-02-15 09:09:55EST ode-example.lisp>
 ;; $Id$
 
 ;;; van der Pol as given in Section 25.5 of the GSL manual.  To
@@ -41,7 +41,7 @@
       (loop
 	 (when (or (>= time max-time) (> iter *max-iter*))
 	   (return (values iter time dep0 dep1)))
-	 make-next-step
+	 (next-step)
 	 (incf iter)
 	 (when print-steps
 	   (format t "~12,6f~10t~12,6f~24t~12,6f~&" time dep0 dep1))))))
diff --git a/ordinary-differential-equations/ode-system.lisp b/ordinary-differential-equations/ode-system.lisp
index 1bb3072578977b3ea0b3012f9e6985c0406b5b6d..4d1f005db2771985a3367f20be40e76d8e51341a 100644
--- a/ordinary-differential-equations/ode-system.lisp
+++ b/ordinary-differential-equations/ode-system.lisp
@@ -1,6 +1,6 @@
 ;; ODE system setup
 ;; Liam Healy, Sun Apr 15 2007 - 14:19
-;; Time-stamp: <2009-02-15 08:42:48EST ode-system.lisp>
+;; Time-stamp: <2009-02-15 09:09:55EST ode-system.lisp>
 ;; $Id$
 
 (in-package :gsl)
@@ -36,8 +36,8 @@
 	    (,step-size (maref ,cstep 0))
 	    ,@(loop for symb in dependent
 		 for i from 0
-		 collect `(,symb (maref ,dep ,i)))
-	    (make-next-step
-	     (apply-evolution
-	      evolve control stepperobj ,ctime ,max-time ,cstep ,dep)))
-	 ,@body))))
+		 collect `(,symb (maref ,dep ,i))))
+	 (flet ((next-step ()
+		  (apply-evolution
+		   evolve ,ctime ,dep ,cstep control stepperobj ,max-time)))
+	   ,@body)))))
diff --git a/ordinary-differential-equations/stepping.lisp b/ordinary-differential-equations/stepping.lisp
index 453b496fff06adfc4d8a95ec289da2a7ef953d83..0ce8ebe3a0f89e2dc4984869b0e69980cf635106 100644
--- a/ordinary-differential-equations/stepping.lisp
+++ b/ordinary-differential-equations/stepping.lisp
@@ -1,6 +1,6 @@
 ;; Stepping functions for ODE systems.
 ;; Liam Healy, Mon Sep 24 2007 - 21:33
-;; Time-stamp: <2009-02-15 08:36:27EST stepping.lisp>
+;; Time-stamp: <2009-02-15 09:19:47EST stepping.lisp>
 ;; $Id$
 
 (in-package :gsl)
@@ -21,6 +21,36 @@
   :initialize-args nil
   :singular (dimension))
 
+#|
+This description applies when scalars=t:
+
+Setup functions for ODE integrators.  The variables `function' and
+`jacobian' are CL functions that define the ODE and its jacobian,
+and should be defined previously with defuns.
+
+The ODE integrator solves a system of first order differential
+equations:
+
+ d y_i / dx = f_i(x,y_1,y_2,...) with y_i=y_1, y_2, ...
+
+The arguments to `function' and `jacobian' aree same: x and y_i as scalars.
+
+Both `function' and `jacobian' return multiple values.  `function'
+returns `dimension' f_i values and `jacobian' returns
+`dimension^2+dimension' values.
+
+`jacobian' returns derivatives of `function' with respect to `x' and y_i.
+The ordering of the `jacobian' values is as follows (the matrix view
+is for convenience only, the values argument is a long list of
+values):
+
+ (values `d f_1 / d x'   `d f_2 / d x'   ... `d f_N / d x'
+         `d f_1 / d y_1' `d f_1 / d y_2' ... `d f_1 / d y_N'
+         `d f_2 / d y_1' `d f_2 / d y_2' ... `d f_2 / d y_N'
+                               ...
+         `d f_N / d y_1' `d f_N / d y_2' ... `d f_N / d y_N')
+|#
+
 (def-make-callbacks ode-stepper
     (function jacobian dimension &optional (scalars t))
   (if scalars
@@ -64,8 +94,8 @@
   "The order of the stepping function on the previous
   step, which can vary if the stepping function itself is adaptive.")
 
-(defmfun step-apply
-    (stepper time step-size y yerr dydt-in dydt-out)
+(defmfun apply-step
+    (stepper time y step-size yerr dydt-in dydt-out)
   "gsl_odeiv_step_apply"
   (((mpointer stepper) :pointer)
    (time :double)
@@ -78,15 +108,15 @@
   :documentation			; FDL
   "Apply the stepping function stepper to the system of equations
    defined by make-ode-stepper, using the step size step-size to
-   advance the system from time time and state y to time t+h.  The new
-   state of the system is stored in y on output, with an estimate of
-   the absolute error in each component stored in yerr If the argument
-   dydt-in is not null it should point an array containing the
-   derivatives for the system at time t on input. This is optional as
-   the derivatives will be computed internally if they are not
-   provided, but allows the reuse of existing derivative information.
-   On output the new derivatives of the system at time t+h will be
-   stored in dydt-out if it is not null.
+   advance the system from time time and state y to time t +
+   step-size.  The new state of the system is stored in y on output,
+   with an estimate of the absolute error in each component stored in
+   yerr If the argument dydt-in is not null it should point an array
+   containing the derivatives for the system at time t on input. This
+   is optional as the derivatives will be computed internally if they
+   are not provided, but allows the reuse of existing derivative
+   information.  On output the new derivatives of the system at time
+   t + step-size will be stored in dydt-out if it is not null.
 
    User-supplied functions defined in the system dydt
    should signal an error or return the correct value.")