Constraining cosmological simulations with peculiar velocities: a forward-modeling approach

TL;DR

The paper addresses how to bridge the gap between cosmological simulations and the uniquely detailed Local Universe by constraining initial conditions from peculiar-velocity data. It introduces Hamlet-PM, a fully Bayesian forward-modeling approach that uses a differentiable particle-mesh gravity solver to infer high-redshift density modes from z=0 velocity measurements, generating 100 constrained DM-only initial conditions evolved with AREPO. The method improves the realism of density and velocity fields in the constrained volume and enables direct comparison to Local-Universe observables, including a quantitative framework (LUM and opt-LUM) for cluster recovery. However, it also reveals tensions, such as an over-abundance of massive halos in the constrained region and an overestimation of velocity amplitudes, likely due to solver limitations, highlighting the need for higher-resolution gravity models and improved observational systematics. Overall, Hamlet-PM provides a flexible, Bayesian platform for constraining cosmological simulations from velocity data and can be extended to incorporate better physics and data in future work.

Abstract

Numerical simulations are a key tool to decipher the dynamics of gravitation. Yet, they fail to spatially reproduce the Universe we observe, limiting comparison between observations and simulations to a statistical level. This is highly problematic for rare, faint or well studied nearby objects that are observed in a single environment. The computational cost of recovering this environment in random simulations is prohibitive. We present Hamlet-PM, a method that enables the constraining of initial conditions for cosmological simulations so as to produce evolved numerical universes that can be directly compared to observations of the Local Universe: constrained simulations. Our method implements the field-level forward modeling of the early-time density field from sparse and noisy measurements of late-time peculiar velocities. The dynamics are integrated with a particle-mesh gravity solver, thus probing the mildly non-linear regime. The code is applied to the Cosmicflows-4 compilation of peculiar velocities up to z < 0.05 (160 Mpc/h). The constrained ICs a re-simulated with a high precision N-body code. A series of one hundred dark-matter only cosmological constrained simulations with a resolution of 512^3 particles in a 500^3 [Mpc/h]3 box is presented. Special attention is given to twelve prominent nearby galaxy clusters, whose simulated counterparts are matched on criteria of mass and separation. We provide a mass estimate constrained by the dynamical environment for each cluster. Field-level forward modeling of the initial conditions produces highly constrained cosmological simulations. Currently, this method already overtakes in quality the pipeline in use in the peculiar-velocity community, although systematic biases still need to be addressed. Furthermore, improving the model is easy thanks to the inherent flexibility of the Bayesian approach.

Figures (12)

• Figure 1: Slices of the mean density field through the super-galactic planes SGX=0 and SGZ=0, cropped to the constrained region. The mean density is computed as a geometrical mean, i.e.$\tilde{\rho} = \mathopen{}\mathclose{\left(\prod_i \rho_i\right)^{1/n}}$. This approach increases the readibility and enables comparison with Hoffman2018Pfeifer2023.
• Figure 2: Power spectrum of the initial as well as the evolved density fields compared with the linear power spectrum. For the reconstructed fields, the solid lines and the shades respectively represent the median and $68\%$ interval around it over the 100 constrained simulations. The initial (resp. evolved) density fields have a resolution of $128^3$ (resp.$512^3$) nodes.
• Figure 3: Halo mass function of the entire volume, in the constrained volume and in the remaining, unconstrained volume. The shades indicate the 1$\sigma$ confidence level derived from the 100 constrained simulations. The empirical halo mass function of Tinker2008 is shown twice; once with the cosmology of the simulation suite and once with the standard Planck 15 values PlanckCollaboration2016. The lower panel shows the ratio to Tinker2008 with the cosmology of the simulations suite.
• Figure 4: Detection rate of major (or nearby) clusters of the Local Universe, defined in Pfeifer2023 as the rate of simulations meeting the upper p-value threshold for each given cluster. Results for the LUM approach are given in blue and results for the opt-LUM algorithm (see \ref{['eq:optlum']}) in red.
• Figure 5: Posterior separation of the matched halos to the observational cluster. The left, central and right columns respectively represent the 3D, radial and ortho-radial separations. The blue and red curves respectively correspond to halos matched with LUM and opt-LUM methods. Values in green indicate the separation between the halos matched by the opt-LUM method and the associated shifted target. The ortho-radial component is normalized by a factor $\sqrt{2}$ so as to enable comparison with the radial component (for an isotropic scatter, $|s_t| / |s_r| = \sqrt{2}$). Clusters are sorted by increasing distance. The separations presented in McAlpine2025, all below $1\,h^{-1}\,{\rm Mpc}$, are not plotted here.