The Arepo public code release

Rainer Weinberger; Volker Springel; Rüdiger Pakmor

Source-linked AI summary

The Arepo public code release

Rainer Weinberger, Volker Springel, Rüdiger Pakmor

arXiv:1909.04667v2Sep 10, 2019astro-ph.IMastro-ph.COastro-ph.GAphysics.comp-ph

TL;DR

Astrophysical simulations require scalable, sophisticated, and well-organized software, while direct collisionless-force calculations have unfavorable N^2 cost. Arepo addresses these needs with a public, documented moving-mesh MHD code coupled to tree-particle-mesh gravity, local time integration, and supporting physics and analysis tools. The release is intended to support reproducibility and further computational astrophysics, although its proposed community development model remains to be tested.

Problem
Astrophysical simulation software must keep pace with hardware and numerical advances, organize complex collaborative development, and avoid the N^2 cost of direct collisionless-force calculations.
Method
Arepo uses finite-volume MHD on a moving unstructured Voronoi mesh, tree-particle-mesh gravity, local time steps, and documented modular code with example physics and analysis tools.
Results
The public release provides a documented, extensible codebase with local time integration, gravitational and MHD solvers, example sub-resolution physics, tests, and usability support.
Takeaways & Limitations
The release is intended to support reproducibility, scientific progress, and community development in computational astrophysics.
Takeaways & Limitations
The MHD 8-wave formulation has known shortcomings, including higher numerical reconnection rates than Athena++ in the magnetic current-sheet test.

Abstract

from arXiv · show

We introduce the public version of the cosmological magnetohydrodynamical moving-mesh simulation code Arepo. This version contains a finite-volume magnetohydrodynamics algorithm on an unstructured, dynamic Voronoi tessellation coupled to a tree-particle-mesh algorithm for the Poisson equation either on a Newtonian or cosmologically expanding spacetime. Time-integration is performed adopting local timestep constraints for each cell individually, solving the fluxes only across active interfaces, and calculating gravitational forces only between active particles, using an operator-splitting approach. This allows simulations with high dynamic range to be performed efficiently. Arepo is a massively distributed-memory parallel code, using the Message Passing Interface (MPI) communication standard and employing a dynamical work-load and memory balancing scheme to allow optimal use of multi-node parallel computers. The employed parallelization algorithms of Arepo are deterministic and produce binary-identical results when re-run on the same machine and with the same number of MPI ranks. A simple primordial cooling and star formation model is included as an example of sub-resolution models commonly used in simulations of galaxy formation. Arepo also contains a suite of computationally inexpensive test problems, ranging from idealized tests for automated code verification to scaled-down versions of cosmological galaxy formation simulations, and is extensively documented in order to assist adoption of the code by new scientific users.

1. INTRODUCTION

Arepo is publicly released to address the need for accessible, documented, verifiable astrophysical simulation software. The release simplifies and documents the code so scientists can use and extend it without direct author support.

Astrophysical simulation progress depends on software and numerical-method advances alongside improvements in computer hardware.
Increasing software complexity creates a need for documentation, modularity, portability, and continuous verification across collaborative teams.
Arepo is a multipurpose gravitational and finite-volume MHD code for cosmic structure formation and broader astrophysical problems.
The code had previously been available to some external users but not to the full astronomical community.
The public release provides a simplified, completely documented version intended for independent use and extension by scientists.It also serves as an overview and introduction for new users, linking equation variables to code parameter names where appropriate.

2. EQUATIONS

Arepo formulates collisionless dynamics, hydrodynamics, and MHD in either static or uniformly expanding cosmological coordinates. Its equations use co-moving variables and include a Friedmann-based expansion model with gravity and source terms.

Arepo solves collisionless-particle and hydrodynamic equations on a uniformly expanding, flat FLRW spacetime.Static spacetime is recovered by setting a = 1 and its time derivatives to zero.
Co-moving variables absorb some effects of expansion, while transformations also modify time and spatial derivatives.
The scale factor evolves through the Friedmann equation under a Λ cold dark matter cosmology, neglecting radiation at relevant redshifts.
The expanding-space gravitational potential is obtained from a co-moving Poisson equation connected to the Friedmann equation and the mean density.
Collisionless phase-space evolution is represented by discrete particles whose trajectories follow equations of motion derived using characteristics.
The ideal-MHD formulation reduces source terms by using time derivatives of auxiliary variables, while external heating and cooling enter the energy equation.

3. MODELING GRAVITATIONAL INTERACTIONS VIA N-BODY DYNAMICS

Arepo combines hierarchical tree and grid-based particle-mesh methods to compute gravitational forces without the N^2 cost of direct particle interactions. Its tree algorithm approximates distant groups while controlling errors and supports softened, periodic, and hybrid force calculations.

Direct collisionless force calculation scales as N^2, so Arepo combines an oct-tree method with a Cartesian-grid Fourier method.
Hierarchical multipole expansion with an oct-tree algorithm: The oct-tree recursively partitions space, computes node multipole moments, and evaluates accepted nodes instead of every individual particle.
Hierarchical multipole expansion with an oct-tree algorithm: Arepo’s relative opening criterion uses node mass, distance, prior acceleration, and α to limit individual interaction errors relative to the total force.
Hierarchical multipole expansion with an oct-tree algorithm: A geometric opening criterion is more costly at equal force accuracy, so the relative criterion is preferred after the first force calculation.
Tree force errors can accumulate into momentum nonconservation, which Arepo mitigates by randomizing tree-domain placement or shifting the periodic box.
Gravitational softening: Gravitational softening suppresses artificial close-particle interactions, with gas-cell softening tied to Voronoi-cell volume through an input parameter.Interactions with different softenings are symmetrized using the larger softening length.
Tree-particle-mesh approach: Tree-particle-mesh calculations discard distant tree interactions beyond a cutoff because the long-range force is supplied by the particle-mesh algorithm.The split scale is specified relative to the PM mesh cell size, and the approach is described as fast, accurate, and flexible.

4. DISCRETIZATION OF MHD

Arepo discretizes ideal MHD with a second-order finite-volume method on a dynamic Voronoi mesh, using reconstructed interface states and fluxes. It combines least-square gradients, divergence cleaning, softened cell-centered gravity, and operator-split gravitational source terms.

MHD discretization: Arepo uses a second-order accurate finite-volume discretization with volume-averaged primitive variables stored at cell centers.Gradients from neighboring cells support piecewise-linear reconstruction, and reconstructed primitive variables are extrapolated to mesh interfaces for flux calculation.
Gradient estimate: The improved gradient estimate minimizes weighted residuals between extrapolated values and neighboring-cell values.The weights are proportional to interface area divided by the distance between cell centers.
Divergence constraint: Powell divergence cleaning advects numerically induced magnetic-field divergences away while remaining flexible for arbitrary mesh geometries and local time steps.The cell divergence uses face magnetic fields, face normals, face areas, and the cell volume.
Divergence constraint: The adopted 8-wave divergence-cleaning formulation has known shortcomings for inclined 2D MHD shock tubes and higher numerical reconnection rates than state-of-the-art Cartesian-grid codes.These limitations are reported alongside its computational flexibility and suitability for cosmological galaxy-formation simulations.
Gravitational interactions of the fluid: Arepo computes gas gravity with cell-centered softened forces whose softening lengths are tied to cell radii.This makes the gravitational field of the gas approximate a space-filling continuous density distribution.
Gravitational interactions of the fluid: Gravity is coupled to the MHD equations by applying half a gravitational source step before and after a full conservative fluid step.For self-gravitating fluids, this operator splitting does not manifestly conserve total energy.

5. COMPUTATIONAL MESH

Arepo uses a fully adaptive, dynamic Voronoi mesh whose generating points move with the fluid and can be locally refined or de-refined. Mesh regularization, distributed-domain ghost points, boundary treatments, and mass-based refinement support spatial adaptivity and robust mesh evolution.

Dynamic Voronoi mesh: Each Voronoi cell is defined by the nearest mesh-generating point, producing interfaces normal to connecting points and equidistant from them.Neighbor connections form the Delaunay tessellation used to construct the mesh.
Parallel mesh construction: In distributed-memory runs, imported ghost points from neighboring domains are required to construct correct Voronoi cells at domain boundaries.Cell ownership follows the generating-point location rather than the cell boundary.
Mesh motion: Mesh-generating points normally follow cell bulk velocities while incorporating pressure-gradient, Lorentz-force, and gravitational accelerations, yielding quasi-Lagrangian behavior.This motion approximately maintains equal-mass cells but does not by itself guarantee regular geometry.
Mesh regularization: Velocity corrections steer highly distorted cells toward their centers of mass when a face-angle distortion threshold is exceeded.The correction uses a fraction of the cell’s characteristic speed, with β controlling the threshold parameterization.
Refinement and de-refinement: Arepo supports local refinement and de-refinement, splitting or removing cells while conservatively redistributing their conserved quantities.Refinement divides quantities by volume ratio, whereas de-refinement distributes them according to neighboring volume overlap.
Refinement and de-refinement: The public version’s mass criterion refines cells above twice the target mass and de-refines cells below half the target.This keeps mass resolution within a narrow corridor around the target, while distorted cells are normally excluded from refinement.
Boundary conditions: Arepo supports periodic, reflective, and inflow/outflow hydrodynamic boundaries, plus periodic or nonperiodic gravitational boundaries.Reflective boundaries mirror mesh-generating points and fluid states so corresponding Riemann solvers return zero mass flux.

6. ADDITIONAL PHYSICS

Arepo includes simple illustrative modules for additional physics beyond gravity and MHD: primordial radiative cooling and unresolved star formation. These examples are documented as guidelines for implementing more sophisticated treatments.

Additional physics: The additional-physics modules are intentionally simple because their system-dependent implementation and verification are less clear than for gravity and MHD.They are included primarily as illustrative guidelines for implementing more sophisticated modules.
Additional physics: The public release includes simple radiative-cooling and star-formation modules as examples of additional physics.The authors note that both models have been widely used in galaxy-formation simulations, although more sophisticated treatments are common in recent work.
Radiative cooling: Radiative cooling is implemented as an energy sink using first-order implicit integration to avoid the tight time-step constraints of explicit integration.The implicit update solves for the final specific internal energy with an iterative root-finding procedure.
Radiative cooling: The public release includes only a primordial H-and-He cooling network under collisionless ionization equilibrium, while extensions are described as straightforward.The cooling rate depends strongly on metallicity and density, making the equation stiff to integrate numerically.
Modeling star formation: The star-formation example follows the Springel-Hernquist model, representing unresolved interstellar gas as cold clouds embedded in a hot, pressure-equilibrated medium.Star-forming cells use a smooth, spatially averaged model rather than resolving the multiphase interstellar-medium structure.
Modeling star formation: Star formation is restricted by temperature and density thresholds, with star particles spawned stochastically according to the cell’s star-formation rate.If insufficient gas remains after an event, the cell is converted fully into a star particle; otherwise, its mass, energy, and momentum are reduced.

7. TIME INTEGRATION

Arepo integrates gravity, MHD, and source terms with separate local constraints organized in nested power-of-two time bins. Its hydrodynamic scheme uses two flux evaluations on a moving mesh, while gravity uses leapfrog integration and active-particle calculations.

Local time stepping: Each cell applies the most restrictive local time-step constraint among hydrodynamics, gravity, and source or sink terms.Dark-matter particles have individual gravitational constraints, while gas cells can also receive MHD, star-formation, and nonlocal wave-arrival constraints.
Local time stepping: Gas-cell hydrodynamic time steps use a CFL criterion whose moving-mesh signal speed accounts for sound and Alfvén speeds.Static meshes additionally include the gas bulk velocity, producing a more restrictive criterion than the moving-mesh case.
Local time stepping: Arepo organizes local steps in a power-of-two hierarchy, allowing elements to use the largest sufficiently small step while keeping shorter steps nested within longer ones.Transitions to longer time steps occur only when the current time is synchronized with the target bin.
Local time stepping: Figure 3 shows each time bin undergoing primitive-variable updates, gradient calculations, two flux calculations, and mesh determination at the step endpoint.Interface fluxes between differently timed cells update both conserved states during the smaller step, leaving larger-bin cells partially updated.
Gravitational integration: Gravity uses a second-order leapfrog scheme alternating position-changing drift operations with velocity-changing kick operations.For fixed step sizes, this produces a symplectic integration scheme.
Gravitational integration: Hierarchical gravity addresses momentum-conservation and performance problems that arise when active particles interact with the full mass distribution under local time stepping.The alternative avoids requiring a full tree construction when only a small fraction of particles needs force calculations.
Hydrodynamic integration: The moving-mesh hydrodynamic update combines beginning-of-step and end-of-step fluxes, with the Riemann problem evaluated twice per time step.The scheme uses time-extrapolated quantities and accounts for time-dependent moving-mesh facet areas.

8. DOMAIN DECOMPOSITION AND LOAD BALANCING

Arepo uses MPI domain decomposition with dynamic balancing of computational and memory loads. It oversamples MPI ranks with multiple spatial chunks and assigns cost factors that reflect hydrodynamic, gravitational, activity, and memory demands.

Load balancing: Arepo uses MPI and dynamically balances computational and memory loads across distributed-memory nodes and processes.This distribution is necessary for large computations on modern supercomputers.
Domain decomposition: The domain decomposition is repeated when the active-particle fraction reaches a specified threshold, balancing gravity, hydrodynamics across time bins, and memory load.The decomposition can be refreshed every time step when enough particles are active.
Load balancing: Oversampling MPI ranks with multiple domain chunks provides additional mappings that smooth residual imbalances across computational work and memory.A rank can combine a computationally expensive, particle-light chunk with a particle-heavy, lower-work chunk.
Domain decomposition: Arepo refines top-level oct-tree leaves whose load exceeds a target fraction, then combines adjacent leaves into similarly loaded domain chunks.The refinement factor is usually set to approximately 2.5-4.0 to create sufficiently fine load discretization.
Load balancing: The secondary balancing step enabled by multiple chunks per MPI rank was not implemented in Gadget-2.
Load estimation: Hydrodynamic and gravitational costs are estimated with separate particle and cell cost factors, including activity frequency and gravity interaction counts.Gravity costs vary with clustering, whereas hydrodynamic cells have a similar average number of neighbors.

9. OTHER CODE FEATURES

Arepo provides practical tools for identifying structures, converting initial conditions, handling large simulation outputs, restarting runs, managing memory, and maintaining portability. These features support long-running, distributed simulations while preserving reproducibility and operational reliability.

Structure and substructure finding: Friends-of-friends and Subfind algorithms identify groups and subhalos either during a run or in post-processing.Subhalos are retained after gravitational unbinding when a bound particle set exceeds the minimum size, and nested subhalos are supported.
Initial conditions conversion: Arepo converts SPH initial conditions into grid-based conditions using an oct-tree tessellation and conservative distribution of SPH quantities.The resulting mesh avoids the highly irregular initial mesh produced by directly using SPH particles.
Data input and output: Snapshot data support HDF5 and older Gadget-compatible formats, with inputs and outputs distributed over arbitrary file counts subject to output limits.Arepo uses explicit POSIX file writes rather than parallel HDF5 or MPI-IO to reduce library-related memory side effects.
Restarting simulations: Checkpointing writes per-rank restart files that reproduce binary-identical results after resumption, provided the same number of MPI ranks is retained.Restart files are written automatically at specified times and in regular intervals, with two recent copies retained.
Memory and implementation: Arepo emphasizes explicit memory control and low fragmentation because cosmological simulations are often limited by physical memory rather than CPU time.Its memory manager supports stack-like freeing and movable blocks for cases where they simplify or accelerate the implementation.
Memory and implementation: The code limits external libraries, uses C11 and MPI-2, and targets compiler portability while warning that aggressive commercial-compiler optimization can produce incorrect code.These choices prioritize maintainability and reliability alongside performance across different systems.

10. EXAMPLES

Arepo is evaluated with automated regression tests spanning smooth flows, discontinuities, interacting shocks, MHD Riemann problems, and cosmological examples. The tests generally recover reference or exact solutions, while the MHD shocktube requires substantially higher resolution than Athena++.

Verification suite: Arepo includes inexpensive example simulations that serve as automated regression tests for functionality and correctness after code changes.The default test set is summarized in Table 7.
Hydrodynamic tests: The sound-wave test probes smooth-flow accuracy and convergence by comparing the evolved perturbation after one domain-crossing time with the initial conditions.The setup uses a perturbation amplitude δ = 10^-6 and CCFL = 0.3.
Hydrodynamic tests: The shock-tube test examines shock, rarefaction, and contact discontinuities while monitoring total variation to constrain numerical diffusivity.The moving mesh produces a shock oscillation, but the flux limiters keep total variation increases below 1%.
Hydrodynamic tests: The interacting-blast-wave test compares a moving-mesh solution at t = 0.038 with a very high-resolution static-grid reference.Initially equally spaced cells move with the flow, without refinement or de-refinement.
MHD tests: The MHD shocktube solution at t = 0.4 is verified against an exact solution and compared directly with Athena++ using 300 cells.Both codes show similar behavior, but Arepo needs roughly three times higher resolution for comparable accuracy.

10.5. Static polytrope

The static-polytrope and multidimensional shock tests probe hydrostatic balance and strong-shock behavior, while the current-sheet test compares numerical reconnection with Athena++.

10.5. Static polytrope: The static-polytrope test initializes a known n = 1 hydrostatic solution with reflective inner and open outer boundaries.The setup spans r = 0.1 to r = 1 and includes both enclosed gas gravity and a static core.
10.5. Static polytrope: Verification compares the polytrope at t = 1 with its initial conditions while keeping the mesh static.The initial velocity is zero and CCFL = 0.3.
10.6. Magnetic current-sheet problem: The current-sheet problem tests Arepo’s MHD scheme, particularly its numerical reconnection rate, in a two-dimensional periodic domain.Cells move, but refinement and de-refinement are disabled.
10.6. Magnetic current-sheet problem: Arepo exhibits a higher numerical reconnection rate than Athena++ in the current-sheet comparison through t = 30.The comparison uses total magnetic energy and includes magnetic-energy-density slices at t = 1.
10.7. Noh shock problem: The Noh shock problem tests the Riemann solver in two and three dimensions using an inward constant-velocity flow that produces a central outward-moving shock.The post-shock solution and shock position are checked against analytic results.

10.8. Gresho vortex

Arepo tests moving-mesh preservation of stationary and smooth vortices, comparing evolved density fields with initial conditions and studying convergence. A small cosmological volume test additionally checks halo-mass results against a previous solution.

10.8. Gresho vortex: The Gresho vortex tests a stationary two-dimensional flow in which centrifugal forces are balanced by a pressure gradient.Arepo evolves the moving mesh to t = 3 and assesses convergence using the density field.
10.8. Gresho vortex: The Gresho setup uses a periodic box of side length L = 1 with piecewise density, azimuthal velocity, and pressure profiles.The initial radial velocity is zero and refinement and de-refinement are not used.
10.8. Yee vortex: The smooth Yee vortex provides a related two-dimensional moving-mesh test with periodic boundaries and a smooth radial solution.The density field is compared with the initial conditions after evolution to t = 10 to study convergence.
Cosmological example: A gravity-only cosmological volume test uses a 50h^-1 Mpc box and 32^3 particles as a computationally inexpensive comparison against a previous solution.The test focuses particularly on the masses of emerging halos.

10.11. Cosmological volume with gas cooling and star formation

Arepo supports cosmological simulations with cooling, star formation, zoomed resolution, galaxy mergers, and flexible Voronoi-mesh initial conditions. The examples verify these capabilities while exposing practical setup and modeling boundaries.

10.11. Cosmological volume with gas cooling and star formation: A 7.5 h^-1 Mpc box with 323 collisionless particles provides the higher mass resolution needed to resolve gas cooling and star formation.
10.11. Cosmological volume with gas cooling and star formation: Gas cells are added by splitting each dark matter particle, and stars form at the centers of dark-matter overdensities.The setup initially contains 2 × 323 resolution elements.
10.11. Cosmological volume with gas cooling and star formation: The cosmological volume example compares its redshift-dependent result with a previous run because the correct solution is unknown.Its star formation rate substantially exceeds observations because energetic feedback processes driving galactic outflows are absent.
10.12. Cosmological zoom simulation with gravity only: Zoom simulations identify and refine the Lagrangian region whose particles end in the target object, while Arepo balances the resulting clustered workload.
10.14. Galaxy merger, including gas and star formation: The galaxy-merger example tracks existing and newly formed stars through two pericenter passages and the final merger.The galaxy centers of mass provide the trajectories, and the resulting star formation rates are compared with a reference solution.
10.12. Cosmological zoom simulation with gravity only: The zoom test compares subhalo masses with a precomputed reference run using identical initial conditions.The z = 0 distribution, identified subhalos, and subhalo mass function are shown as diagnostics.
10.15. General remarks on initial condition generation: Arepo’s evolving Voronoi coordinates allow resolution-element geometry and local resolution to adjust flexibly to the problem.In approximately spherical configurations, radially adaptive cell volumes can keep the mass per cell approximately equal.
10.15. General remarks on initial condition generation: Voronoi mesh construction requires robust geometric predicates and exact arithmetic for degeneracies, while cell-averaged variables complicate high-quality initial conditions.Perfect Cartesian starting meshes can therefore spend considerable time in early mesh constructions.

11. CODE DEVELOPMENT ASPECTS

Arepo’s growing complexity and user base required a more formal development model. The public release therefore emphasizes a stable, documented base, modular extensions, and community-based support.

11. CODE DEVELOPMENT ASPECTS: Arepo’s nearly 10^5-line public code is about five times larger than Gadget-2 because mesh construction, hydrodynamics, and other numerical components are more sophisticated.
11. CODE DEVELOPMENT ASPECTS: A master branch contains tested modules, feature branches isolate development, and senior-developer approval is required for integration.
11. CODE DEVELOPMENT ASPECTS: The public version uses a stable base intended mainly for bug fixes and examples, with new features developed in branches, forks, or separate repositories.Patches and stable model branches can remain independent of the original authors’ involvement.
11. CODE DEVELOPMENT ASPECTS: Arepo’s resemblance to Gadget-2 and Gadget-3 eases entry for experienced users and simplifies porting existing modules.Module patches should minimize changes to the original code and place added functionality in new source files.
11. CODE DEVELOPMENT ASPECTS: The release prioritizes complete documentation, usage guidance, and small working examples to train users as simulation software becomes more complex.A support forum extends this documentation by enabling users to report and answer problems publicly.

12. SUMMARY AND CONCLUSIONS

Arepo is released publicly as a documented, open-source moving-mesh code with gravity, MHD, local time integration, and supporting analysis and initialization tools. The release provides a foundation for community use and development while identifying important directions for future improvement.

Core algorithms: Arepo combines finite-volume MHD on a moving unstructured Voronoi mesh with tree-particle-mesh gravity for multiple collisionless particle species.This supports simulations of cosmic structure formation and other astrophysical problems.
Time integration and physics: Individual local time-step criteria enable efficient calculations across large dynamic ranges in time.Additional gas physics, including radiative cooling and a simple unresolved interstellar-medium model, can be coupled to the main code.
Supporting tools: Built-in FOF, Subfind, and SPH-to-Arepo conversion routines support productive use, including on-the-fly or post-processing analysis.These tools are especially aimed at simulations of cosmic structure formation.
Community development: The release provides completely documented source code, examples, and a support forum to facilitate adoption and community-based extension.Developers may contribute modules through independent repositories, patches, or public branches.
Future development: Further work is needed to improve higher-order methods, additional physical modules, distributed-memory scalability, GPU support, and energy efficiency.The authors specifically identify aggressive zoom simulations and future supercomputer use as challenging settings.
Scientific role: The authors frame the public release as supporting reproducibility and scientific progress as computational astrophysics software becomes more complex.They present continuous improvement of numerical methods and codes as a prerequisite for turning hardware advances into scientific answers.