A Stable, Accurate Methodology for High Mach Number, Strong Magnetic Field MHD Turbulence with Adaptive Mesh Refinement: Resolution and Refinement Studies

Performing a stable, long duration simulation of driven MHD turbulence with a high thermal Mach number and a strong initial magnetic field is a challenge to high-order Godunov ideal MHD schemes because of the difficulty in guaranteeing positivity of the density and pressure. We have implemented a robust combination of reconstruction schemes, Riemann solvers, limiters, and Constrained Transport EMF averaging schemes that can meet this challenge, and using this strategy, we have developed a new Adaptive Mesh Refinement (AMR) MHD module of the ORION2 code. We investigate the effects of AMR on several statistical properties of a turbulent ideal MHD system with a thermal Mach number of 10 and a plasma $\beta_0$ of 0.1 as initial conditions; our code is shown to be stable for simulations with higher Mach numbers ($M_rms = 17.3$) and smaller plasma beta ($\beta_0 = 0.0067$) as well. Our results show that the quality of the turbulence simulation is generally related to the volume-averaged refinement. Our AMR simulations show that the turbulent dissipation coefficient for supersonic MHD turbulence is about 0.5, in agreement with unigrid simulations.

💡 Research Summary

This paper presents a robust numerical framework for conducting long‑duration simulations of driven, supersonic magnetohydrodynamic (MHD) turbulence under conditions typical of star‑forming molecular clouds: high thermal Mach numbers (up to M ≈ 17) and low plasma beta (β ≈ 0.1 down to 0.0067). The authors identify the principal obstacle in such regimes as the inability of high‑order Godunov schemes to guarantee positivity of density and pressure, which leads to numerical blow‑up when strong shocks and strong magnetic fields coexist. To overcome this, they assemble a carefully chosen combination of reconstruction methods (PLM or PPM), slope limiters (from highly diffusive minmod to less diffusive monotized central), Riemann solvers (Roe, HLL, HLLD), and, most critically, a Constrained Transport (CT) scheme with three alternative electromotive‑force (EMF) averaging strategies (simple arithmetic, face‑to‑edge integration, and flux‑sign‑based selection). This CT‑based approach enforces ∇·B = 0 to machine precision and, unlike divergence‑cleaning methods, does not compromise pressure positivity.

The implementation is built into the ORION2 code, which uses the Chombo block‑structured adaptive mesh refinement (AMR) library. The AMR hierarchy follows the Berger‑Colella algorithm, refined by a factor of two at each level, with up to three refinement levels (effective resolution 1024³). Time sub‑cycling is employed so that each level advances with a timestep reduced by the same factor as its spatial refinement. Synchronization steps include flux correction for conserved variables and a dedicated magnetic‑field correction that ensures the same electric field is used on coarse‑fine interfaces, preserving the divergence‑free condition across levels. Interpolation of face‑centered magnetic fields uses a two‑step process (coarse‑to‑fine face interpolation followed by a projection step based on Martin & Colella 2000) that maintains discrete solenoidality without requiring the ghost fields themselves to be divergence‑free.

A suite of tests validates the method. The primary scientific runs adopt an initial uniform magnetic field corresponding to β = 0.1 and drive turbulence to a statistically steady state with rms Mach number ≈ 10. A more extreme case with M ≈ 17.3 and β ≈ 0.0067 demonstrates the scheme’s stability at the limits of current computational capability. In all cases the density and pressure remain strictly positive, and the turbulent dissipation coefficient ε (defined via the decay of kinetic energy) settles around 0.5, matching results from uniform‑grid simulations at the same effective resolution.

A central focus of the study is the relationship between the fraction of the domain that is refined (the “refinement coverage”) and the fidelity of turbulence statistics. By varying refinement criteria (based on density, velocity gradients, and magnetic‑field gradients) the authors produce runs with refinement coverages ranging from ≈ 30 % to ≈ 90 % of the volume at the finest level. They find that when the volume‑averaged refinement exceeds roughly 70 %, the velocity power spectrum retains the expected –2 to –2.2 slope, the probability‑density functions (PDFs) of density and magnetic field exhibit the same high‑density and low‑density tails as the corresponding uniform‑grid runs, and the measured ε remains unchanged. Below this threshold, the dissipation rate is underestimated, the power‑law slope steepens artificially, and the PDFs become biased, indicating that insufficient refinement degrades the representation of small‑scale structures and shock fronts.

Performance-wise, the AMR runs achieve comparable statistical accuracy to a 1024³ uniform grid while using only about 40 % of the total cells, translating into substantial savings in CPU time and memory. The authors also discuss the overhead associated with AMR synchronization and magnetic‑field correction, noting that it remains modest relative to the gains from reduced resolution in low‑density regions.

In conclusion, the paper demonstrates that a carefully engineered combination of high‑order Godunov reconstruction, appropriate Riemann solvers, robust limiters, and a divergence‑preserving CT scheme enables stable, accurate simulations of highly supersonic, low‑beta MHD turbulence on adaptive meshes. This capability opens the door to realistic star‑formation simulations that can resolve collapsing cores with high fidelity while efficiently treating the surrounding diffuse medium. Future work will integrate self‑gravity, radiative transfer, and more sophisticated refinement criteria to further broaden the applicability of the method.