Second-order accurate finite volume methods are among the most popular algorithms for the simulation of astrophysical phenomena that can be approximated by the ideal magnetohydrodynamics (MHD) equations. While such methods are robust and offer decent accuracy with relative simplicity, higher-order methods with increased arithmetic intensity may be more efficient on modern architectures with poor machine balance (ratio of floating-point operations per cycle to sustained memory operations per cycle).
We present the results of a fourth-order accurate finite volume method for ideal MHD developed for the Athena++ astrophysics code. The method employs a constrained transport (CT) formulation to prevent the onset of numerical monopoles. Challenging multidimensional astrophysics problems are used as benchmarks to compare the accuracy, robustness and computational efficiency of the fourth-order method to the existing second-order scheme used in Athena++. We examine the arithmetic intensity and data locality tradeoffs in vectorized, parallel (MPI) domain-decomposed simulations performed on several multicore architectures, including Intel Xeon Phi Knights Landing processors.