def rk5_stepper(x, h, y , deriv, errors=False):
    """ Compute y(x+h) from y and dy/dx=deriv(x,y).

    Uses a one-step 5th-order Runge-Kutta algorithm.

    Returns x+h, y(x+h) if errors is False; otherwise
    returns x+h, y(x+h), yerr where yerr is an error
    estimate.

    Adapted from Numerical Recipes.
    """
    k1 = h * deriv(x,y)
    k2 = h * deriv(x+0.2*h, y+0.2*k1)
    k3 = h * deriv(x+0.3*h, y+0.075*k1+0.225*k2)
    k4 = h * deriv(x+0.6*h, y+0.3*k1-0.9*k2+1.2*k3)
    k5 = h * deriv(x+h, y-.2037037037037037037037037*k1
             +2.5*k2-2.592592592592592592592593*k3
             +1.296296296296296296296296*k4)
    k6 = h * deriv(x+0.875*h, y+.2949580439814814814814815e-1*k1
             +.341796875*k2+.4159432870370370370370370e-1*k3
             +.4003454137731481481481481*k4+.61767578125e-1*k5)
    yn = y + (.9788359788359788359788361e-1*k1
             +.4025764895330112721417070*k3
             +.2104377104377104377104378*k4
             +.2891022021456804065499718*k6)
    xn = x+h
    if errors:
        yerr = (-.429377480158730158730159e-2*k1
                +.186685860938578329882678e-1*k3
                -.341550268308080808080807e-1*k4
                -.1932198660714285714285714e-1*k5
                +.391022021456804065499718e-1*k6)
        return xn,yn,yerr
    else:
        return xn,yn
Compute y(x+h) from y and dy/dx=deriv(x,y).

    Uses a one-step 5th-order Runge-Kutta algorithm.

    Returns x+h, y(x+h) if errors is False; otherwise
    returns x+h, y(x+h), yerr where yerr is an error
    estimate.

    Adapted from Numerical Recipes.