An optimal, fully implicit algorithm for the low-$\beta$ ``two-field'' two-fluid MHD model.

The low-$\beta$ ``two-field'' two-fluid magnetohydrodynamics model is appealing owing to its simplicity and its wide applicability in strongly magnetized plasmas. However, it supports fast dispersive waves that challenge its numerical integration, and demand efficient implicit integration methods. We propose an efficient, parallel, optimal nonlinearly implicit algorithm for the low-$\beta$ XMHD model based on physics-based preconditioned Jacobian-free Newton-Krylov (JFNK) methods.\footnote{L. Chac\'on and A. Stanier, \emph{J. Comput. Phys.}, submitted (2016)} The proposed preconditioner leverages earlier developments of effective physics-based preconditioners for MHD\footnote{L. Chac\'on, D. A. Knoll, and J. M. Finn, \emph{J. Comput. Phys.}, \textbf{178} (2002)} and high-$\beta$ extended MHD.\footnote{L. Chac\'on and D. A. Knoll, \emph{J. Comput. Phys.}, \textbf{188} (2003) } We demonstrate the performance of the algorithm with challenging numerical examples. In particular, we demonstrate optimal weak parallel scaling for a fixed implicit time step on grids up to 4096$\times$4096 on 4096 cores (the maximum available to us). At these resolutions, speedups with respect to explicit algorithms reach up to four orders of magnitude.