A solver for linear scalar ordinary differential equations whose running time is bounded independent of frequency

When the eigenvalues of the coefficient matrix for a linear scalar ordinary differential equation are of large magnitude, its solutions exhibit complicated behaviour, such as high-frequency oscillations, rapid growth or rapid decay. The cost of representing such solutions using standard techniques grows with the magnitudes of the eigenvalues. As a consequence, the running times of most solvers for ordinary differential equations also grow with these eigenvalues. However, a large class of scalar ordinary differential equations with slowly-varying coefficients admit slowly-varying phase functions that can be represented at a cost which is bounded independent of the magnitudes of the eigenvalues of the corresponding coefficient matrix. Here, we introduce a numerical algorithm for constructing slowly-varying phase functions which represent the solutions of a linear scalar ordinary differential equation. Our method's running time depends on the complexity of the equation's coefficients, but is bounded independent of the magnitudes of the equation's eigenvalues. Once the phase functions have been constructed, essentially any reasonable initial or boundary value problem for the scalar equation can be easily solved. We present the results of numerical experiments showing that, despite its greater generality, our algorithm is competitive with state-of-the-art methods for solving highly-oscillatory second order differential equations. We also compare our method with Magnus-type exponential integrators and find that our approach is orders of magnitude faster in the high-frequency regime.