Numerical solution of ordinary differential equations in geodetic science using adaptive Gauss numerical integration method

Abstract

In this paper a new method of numerically solving ordinary differential equations is presented. This method is based on the Gaussian numerical integration of different orders. Using two different orders for numerical integration, an adaptive method is derived. Any other numerical solver for ordinary differential equations can be used alongside this method. For instance, Runge–Kutta and Adams–Bashforth–Moulton methods are used together with this new adaptive method. This method is fast, stable, consistent, and suitable for very high accuracies. The accuracy of this method is always higher than the method used alongside with it. Two applications of this method are presented in the field of satellite geodesy. In the first application, for different time periods and sampling rates (increments of time), it is shown that the orbit determined by the new method is—with respect to the Keplerian motion—at least 6 and 25,000,000 times more accurate than, respectively, Runge–Kutta and Adams–Bashforth–Moulton methods of the same degree and absolute tolerance. In the second application, a real orbit propagation problem is discussed for the GRACE satellites. The orbit is propagated by the new numerical solver, using perturbed satellite motion equations up to degree 280. The results are compared with another independent method, the Unscented Kalman Filter. It is shown that the orbit propagated by the numerical solver is approximately 50 times more accurate than the one propagated by the Unscented Kalman Filter approach.