Efficient method for solving highly oscillatory ordinary differential equations with applications to physical systems

Abstract

We present a novel numerical routine (oscode) with a C++ and Python interface for the efficient solution of one-dimensional, second-order, ordinary differential equations with rapidly oscillating solutions. The method is based on a Runge-Kutta-like stepping procedure that makes use of the Wentzel-Kramers-Brillouin (WKB) approximation to skip regions of integration where the characteristic frequency varies slowly. In regions where this is not the case, the method is able to switch to a madeto-measure Runge-Kutta integrator that minimises the total number of function evaluations. We demonstrate the effectiveness of the method with example solutions of the Airy equation and an equation exhibiting a burst of oscillations, discussing the error properties of the method in detail. We then show the method applied to physical systems. First, the one-dimensional, time-independent Schrödinger equation is solved as part of a shooting method to search for the energy eigenvalues for a potential with quartic anharmonicity. Then, the method is used to solve the Mukhanov-Sasaki equation describing the evolution of cosmological perturbations, and the primordial power spectrum of the perturbations is computed in different cosmological scenarios. We compare the performance of our solver in calculating a primordial power spectrum of scalar perturbations to that of BINGO, an efficient code specifically designed for such applications, and find that our method performs better.

Beyond the Runge-Kutta-Wentzel-Kramers-Brillouin method

“…The RKWKB method (to second-order in the WKB series) is obtained by replacing the Taylor series solution with the two independent WKB solutions in (2.6). As the WKB solutions typically remain good approximations to the true solution over many oscillations, unlike truncated Taylor series solutions, the RKWKB method is able to take much larger steps than Runge-Kutta methods, which in general reduces the required run time [5].…”

“…Handley et al [4] and Agocs et al [5] also enable the program to switch to a standard Runge-Kutta-Fehlberg method in slowly oscillating regions, however in this paper we use RKWKB to refer to the simpler non-switching approach.…”

“…This means that, while similar, the Jordan-Magnus approach remains distinct from WKB. In certain senses this lack ofω terms may be an advantage as derivatives can be difficult to evaluate numerically [5].…”

“…Handley, Lasenby and Hobson proposed in [4], and developed with Agocs in [5], a new numerical method based on a combination of the standard Runge-Kutta (RK) approach and the Wentzel-Kramers-Brillouin (WKB) approximation, dubbed the "Runge-Kutta-Wentzel-Kramers-Brillouin" (RKWKB) method. In highly oscillatory regions of the solution using the WKB approximation allows integration steps which span * jarb2@alumni.cam.ac.uk † wh260@mrao.cam.ac.uk many oscillations, giving faster performance than standard RK methods.…”

“…In highly oscillatory regions of the solution using the WKB approximation allows integration steps which span * jarb2@alumni.cam.ac.uk † wh260@mrao.cam.ac.uk many oscillations, giving faster performance than standard RK methods. However RKWKB is limited to essentially one dimensional systems [4,5]. We can write a multidimensional coupled system of linear ODEs in the form of (1.1)…”