Nonlinear differential equation solvers via adaptive Picard–Chebyshev iteration: Applications in astrodynamics

An adaptive self-tuning Picard–Chebyshev numerical integration method is presented for solving initial and boundary value problems by considering high-fidelity perturbed two-body dynamics. The current adaptation technique is self-tuning and adjusts the size of the time interval segments and the number of nodes per segment automatically to achieve near-maximum efficiency. The technique also uses recent insights on local force approximations and adaptive force models that take advantage of the fixed-point nature of the Picard iteration. In addition to developing the adaptive method, an integral quasi-linearization “error feedback” term is introduced that accelerates convergence to a machine precision solution by about a of two. The integral quasi linearization can be implemented for both first- and second-order systems of ordinary differential equations. A discussion is presented regarding the subtle but significant distinction between integral quasi linearization for first-order systems, second-order systems that can be rearranged and integrated in first-order form, and second-order systems that are integrated using a kinematically consistent Picard–Chebyshev iteration in cascade form. The enhanced performance of the current algorithm is demonstrated by solving an important problem in astrodynamics: The perturbed two-body problem for near-Earth orbits. The adaptive algorithm has proven to be more efficient than an eighth-order Gauss–Jackson and a 12th/10th-order Runge–Kutta while maintaining machine precision over several weeks of propagation.