Exponential Integrators for Resistive Magnetohydrodynamics: Matrix-free Leja Interpolation and Efficient Adaptive Time Stepping

We propose a novel algorithm for the temporal integration of the resistive magnetohydrodynamics (MHD) equations. The approach is based on exponential Rosenbrock schemes in combination with Leja interpolation. It naturally preserves Gauss's law for magnetism and is unencumbered by the stability constraints observed for explicit methods. Remarkable progress has been achieved in designing exponential integrators and computing the required matrix functions efficiently. However, employing them in MHD simulations of realistic physical scenarios requires a matrix-free implementation. We show how an efficient algorithm based on Leja interpolation that only uses the right-hand side of the differential equation (i.e., matrix free) can be constructed. We further demonstrate that it outperforms Krylov-based exponential integrators as well as explicit and implicit methods using test models of magnetic reconnection and the Kelvin-Helmholtz instability. Furthermore, an adaptive step-size strategy that gives excellent and predictable performance, particularly in the lenient- to intermediate-tolerance regime that is often of importance in practical applications, is employed.