On the basis of symmetric E-methods with higher derivatives having the convergence order four, six, or eight, implicit extrapolation schemes are constructed for the numerical solution of ordinary differential equations. The combined step size and order control used in these schemes implements an automatic global error control in the extrapolation E-methods, which makes it possible to solve differential problems in automatic mode up to the accuracy specified by the user (without taking into account round-off errors). The theory of adjoint and symmetric methods presented in this paper is an extension of the results that are well known for the conventional Runge-Kutta schemes to methods involving higher derivatives. Since the implicit extrapolation based on multi-stage Runge-Kutta methods can be very time consuming, special emphasis is made on the efficiency of calculations. All the theoretical conclusions of this paper are confirmed by the numerical results obtained for test problems.

CEMAT - Center for Computational and Stochastic Mathematics