MHDISMv4-CT (Magnetohydrodynamic InterStellar Medium version 4 – Constrained Transport) is a Julia-based finite volume code for solving the ideal magnetohydrodynamics (MHD) equations with a divergence-free magnetic field evolution guarantee. The code employs the Constrained Transport (CT) method on a staggered grid (Yee, 1966) to maintain ∇·B = 0 to machine precision. This document provides a comprehensive technical description of the code architecture, numerical algorithms, measured performance benchmarks, parallelisation strategy, and its scientific heritage building upon two decades of interstellar medium (ISM) and Local Bubble simulations.
Scientific Motivation
The interstellar medium (ISM) is a complex, multi-phase environment shaped by stellar feedback, magnetic fields, and turbulence. Accurate numerical modelling of the ISM dynamics requires robust magnetohydrodynamic (MHD) solvers that preserve the divergence-free constraint on the magnetic field (∇·B = 0). Violations of this constraint lead to unphysical forces and numerical instabilities. MHDISMv4-CT represents the latest evolution in a family of MHD codes developed over two decades for studying the ISM, supernova-driven turbulence, and the Local Bubble – the low-density cavity of hot plasma surrounding our Solar System.
Code Background
This code builds upon extensive previous work on ISM simulations, organised into two code families:
- EAF Family (FORTRAN 90/2003, AMR-based)
-
- EAF: Évora Astrophysical Fluid code – A hydrodynamical MPI-based code with Adaptive Mesh Refinement developed by de Avillez for his Ph.D. thesis on Galactic Fountains (de Avillez, 2000) and subsequent research on disk-halo interaction and ISM turbulence (de Avillez & Mac Low, 2001, 2002; de Avillez & Berry, 2001; de Avillez & Breitschwerdt, 2004). EAF featured PPM reconstruction and was designed for large-scale 3D simulations of the supernova-driven ISM.
- EVAF-PAMR: Évora-Vienna Parallel Adaptive Mesh Refinement code – A Fortran 2003 MHD extension of EAF developed by de Avillez while on sabbatical leave at the University of Vienna (2004/2005). EVAF-PAMR added magnetic field evolution using the Dai & Woodward (1994, 1998) solver for ideal MHD, enabling the first 3D MHD ISM-Galactic Fountain scale (de Avillez & Breitschwerdt, 2005a,b, 2007; Mac Low et al., 2005) and Local Bubble (Breitschwerdt & de Avillez, 2006; de Avillez, 2009) simulations.
- EVAF-PAMR NEI Extension: EVAF-PAMR NEI Extension introduced for the first time non-equilibrium ionisation (NEI), time-dependent cooling and emission in disk-halo interaction studies, in the Local Bubble and electron transport in the ISM in collaboration with TU Berlin (de Avillez & Breitschwerdt, 2010, 2012a,b,c; de Avillez et al., 2012, 2020; Breitschwerdt & de Avillez, 2021).
- MHDISM Family Codes (Julia, structured grids)
-
- MHDISM series (v1–v4): A new family of cell-centred MHD codes written in Julia for ISM simulations. Progressive development incorporated ideal, resistive, and Hall MHD physics, multispecies evolution with 112 ions/atoms and electrons, non-equilibrium ionisation, time-dependent cooling, and variable γ parameter—enabling the equation of state to evolve with the ionic structure via the internal energy of the plasma. The multispecies and thermodynamic evolution is based on the Collisional+Photo Ionization Plasma Emission Software (CPIPES; de Avillez 2018; de Avillez et al. 2018; see also de Avillez et al. 2020).
- MHDISMv4.02: The development branch of the cell-centered MHDISM series v1–v4, featuring divergence cleaning methods: GLM (Dedner et al., 2002), EGLM (Dedner et al., 2002), and projection (Brackbill & Barnes, 1980). The code implements primitive, conserved, and characteristic variable reconstruction, along with a variety of reconstruction methods and flux solvers for the MHD equations.
- MHDISMv4-CT: The current code is a new development branch of MHDISMv4.02 implementing the Constrained Transport (CT) method (Evans & Hawley,1988; Balsara,1999; Tóth, 2000) on a staggered grid (Yee, 1966), preserving ∇·B = 0 to machine precision. The scheme uses cell-centered hydrodynamic variables and face-centred magnetic field variables.

