E pur si muove: Galiliean-invariant cosmological hydrodynamical simulations on a moving mesh

TL;DR

This paper introduces AREPO, a moving-mesh hydrodynamics code based on a Galilean-invariant Voronoi tessellation that moves with the flow to provide adaptive resolution without mesh tangling. It couples a second-order Godunov finite-volume solver on a moving unstructured mesh with self-gravity, using improved surface-based gravity work terms to achieve accurate energy conservation and realistic structure formation. Through extensive hydrodynamic and self-gravity test suites, AREPO demonstrates competitive accuracy for shocks, contact discontinuities, and fluid instabilities, while mitigating the Galilean non-invariance issues of Eulerian codes. The results suggest that the moving-mesh approach provides a versatile, scalable alternative to SPH and AMR for cosmological and other complex fluid problems, with broad potential extensions to magnetohydrodynamics and beyond.

Abstract

Hydrodynamic cosmological simulations at present usually employ either the Lagrangian SPH technique, or Eulerian hydrodynamics on a Cartesian mesh with adaptive mesh refinement. Both of these methods have disadvantages that negatively impact their accuracy in certain situations. We here propose a novel scheme which largely eliminates these weaknesses. It is based on a moving unstructured mesh defined by the Voronoi tessellation of a set of discrete points. The mesh is used to solve the hyperbolic conservation laws of ideal hydrodynamics with a finite volume approach, based on a second-order unsplit Godunov scheme with an exact Riemann solver. The mesh-generating points can in principle be moved arbitrarily. If they are chosen to be stationary, the scheme is equivalent to an ordinary Eulerian method with second order accuracy. If they instead move with the velocity of the local flow, one obtains a Lagrangian formulation of hydrodynamics that does not suffer from the mesh distortion limitations inherent in other mesh-based Lagrangian schemes. In this mode, our new method is fully Galilean-invariant, unlike ordinary Eulerian codes, a property that is of significant importance for cosmological simulations. In addition, the new scheme can adjust its spatial resolution automatically and continuously, and hence inherits the principal advantage of SPH for simulations of cosmological structure growth. The high accuracy of Eulerian methods in the treatment of shocks is retained, while the treatment of contact discontinuities improves. We discuss how this approach is implemented in our new parallel code AREPO, both in 2D and 3D. We use a suite of test problems to examine the performance of the new code and argue that it provides an attractive and competitive alternative to current SPH and Eulerian techniques. (abridged)

Figures (50)

• Figure 1: Example of a Voronoi and Delaunay tessellation in 2D, with periodic boundary conditions. The panel on the left shows the Voronoi tessellation for $N=64$ points (shown as red circles), the panel in the middle gives the corresponding Delaunay tessellation, while the panel on the right shows both simultaneously (solid lines show the Voronoi, dashed lines the Delaunay tessellation).
• Figure 2: The point insertion algorithm in 2D. We start with a valid Delaunay triangulation in which we want to insert an additional point. We first locate the triangle containing the point (step 1), then split it into three triangles (step 2). The edges (drawn in red) in the new triangles opposite of the inserted point may violate the in-circle criterion and need to be tested individually. If an edge is Delaunay (step 3), it is part of the final tessellation, but if it violates the in-circle criterion (step 4), the edge needs to be flipped in the quadrilateral formed by the adjacent triangles (step 5). The flip generates additional edges that need to be tested (steps 6 and 7). Any violating edge found (e.g. step 9) needs to be corrected by flips. Once all remaining new edges are validated (steps 10 and 11), we arrive again at a valid Delaunay tessellation (step 12).
• Figure 3: A '1-to-4' flip. A newly inserted point splits its insertion tetrahedron into 4 daughter tetrahedra.
• Figure 4: The standard replacement operation in 3D required for restoring Delaunayhood. It consists of 2-to-3 (from left to right) or 3-to-2 (from right to left) flips of tetrahedra. Note that the 2-to-3 flip is only possible if the line connecting the two points opposite of the common face intersects the interior of this face. Conversely, the 3-to-2 flip is only possible if an edge is shared by exactly three tetrahedra.
• Figure 5: Point insertion in 2D in the normal case (top, via a 1-to-3 flip) and the degenerate case (bottom), where the point lies exactly on an edge of the current tessellation. In the latter case, the two triangles need to be replaced with four triangles (a 2-to-4 flip).