fixed classic rk4 step actually returns result for two half steps

[Krishna Myneni]

The function gsl_odeiv2_driver_apply_fixed_step() function in GSL 2.5,
using the rk4 step method actually returns the result for two half-
steps instead of the result for a single step. The attached code shows
the result for integrating the Lorenz system of ODEs with given
parameters and initial conditions, after 1 step at dt = 1e-4:

$ ./bench_gsl_rk4
1.59864268529822870e-03  9.99903676398808394e-01  7.99388902928489989e-
08

The above result for y[0] only has 12 significant digits in comparison
to the double precision result for a single step of the classical RK4
estimate for the step size used. Of course it is more accurate, but the
point is that it is not giving the expected result for 1 classical RK4
step. The problem appears to be with rk4_apply() which is using both
single step and two half steps for estimating the error. It leaves the
output state to be the result of the two half steps. 

Perhaps this is by design, but it will cause confusion when comparing
with other rk4 fixed step solutions.

Krishna Myneni

bench_gsl_rk4.c
Description: Text Data

signature.asc
Description: This is a digitally signed message part