The Algorithm
The method used for integration is an embedded RK2(3) scheme due to Fehlberg, and reported on
page 170 of Hairer et al. (1987). Lagged variables (and gradients) are stored in a
ringbuffer at each step of the integrator. Interpolation is required to estimate values of the
lagged variables between storage times. For numerical probity it is essential that the
interpolation of lagged variables is of a higher order of approximation than the integrator,
otherwise the assumptions underlying the error estimate from the RK pair will not be met. The
algorithm used in Solv95 uses cubic hermite interpolation (e.g. Burden and Faires 1987) to
achieve this (which is the reason that gradients need to be stored along with lagged values). The
consequences of not using consistent interpolation and integration schemes are vividly
illustrated in Highman (1993). Paul (1992) was also influential in the design of the method used
here, and the step size selection is straight out of Press et al. (1992) (method, not
code!). The RK2(3) pair used is not actually optimal - it should be possible to derive an
improved scheme - see Butcher (1987) for an explanation of how to go about it.
A deficiency of the algorithm is that switches are not tracked automatically - if a switch
impacts on a lagged variable the resulting derivative discontinuities are only dealt with by
adaptive timestepping through them (Note, however, that lagged variables and gradients are stored
immediately before and immediately after the application of a switch - i.e. the interpolator is
not applied blindly through switches!).
- Burden, R.L. and J.D. Faires (1985) Numerical Analysis. Pridle Weber and Schmidt, Boston.
- Butcher, J.C. (1987) The Numerical Analysis of Ordinary Differential
Equations. John Wiley & Sons, Chichester.
- Hairer, E., S.P.Norsett & G.Wanner (1987) Solving Ordinary differential
Equations I. Springer-Verlag Berlin. p170 RKF2(3)B
- Highman, D.J. (1993) Appl. Numer. Math. 12:403-414
- Paul, C.A.H (1992) Appl. Numer. Math. 9:403-414
- Press et al. (1992) Numerical Recipes in C. CUP
Back to index