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

Title: Constraining cosmological simulations with peculiar velocities: a forward-modeling approach (Hamlet-PM)

Problem Statement
Numerical simulations are essential for understanding gravitational dynamics, yet standard cosmological simulations only reproduce the Universe statistically, failing to spatially match the specific "Local Universe" observed by astronomers. This limitation hinders the study of rare, faint, or well-studied nearby objects that exist in a single unique environment. While "constrained simulations" aim to reproduce the local environment, existing pipelines (such as the CLUES consortium's workflow) rely on a heterogeneous, multi-step approach: reconstructing density/velocity fields via a Wiener Filter (WF), shifting constraints back in time using the Reverse Zeldovich Approximation (RZA), and generating initial conditions (ICs). These steps rely heavily on linear theory and assumptions (e.g., Gaussian velocity uncertainties) that break down in the non-linear regime and introduce systematic biases. Furthermore, current methods often struggle to accurately model the dynamics of the late-time Universe where peculiar velocities are measured.

Methodology: Hamlet-PM
The authors present Hamlet-PM, a fully Bayesian, forward-modeling framework that constrains initial conditions directly from late-time ($z=0$) peculiar velocity data. Unlike previous iterations (e.g., Hamlet-L) that used linear theory, Hamlet-PM integrates the following components:

1. Differentiable Gravity Solver: The method replaces the linear model with a coarse particle-mesh (PM) gravity solver (pmwd) implemented in JAX. This allows for the integration of dynamics into the mildly non-linear regime.
2. Hybrid Velocity Field Construction: To address issues with empty cells in void regions where standard Cloud-in-Cell (CiC) interpolation fails, the authors employ a hybrid approach. In dense regions, the velocity field is derived from the CiC representation of the N-body simulation. In nearly empty regions, a linear theory approximation is used. These are blended using a density-dependent weight function.
3. Bayesian Inference: The framework estimates the posterior probability distribution of the high-redshift overdensity field Fourier modes and the distances of the constraints. It utilizes Hamiltonian Monte Carlo (HMC) via the blackjax library to explore the parameter space, avoiding the approximations of the Wiener Filter.
4. Data: The method is applied to the Cosmicflows-4 (CF4) catalog, a composite dataset of $\sim56,000$ galaxies grouped into $\sim38,000$ groups, limited to $cz \leq 16,500$ km/s ($z \sim 0.05$) to ensure data isotropy and manage computational constraints.
5. Simulation Execution: The inferred initial conditions (at low resolution, $128^3$) are augmented with random small-scale power to reach $512^3$ particles. These are evolved to $z=0$ using the AREPO code, generating a suite of 100 constrained simulations in a $500^3$ $h^{-1}\text{Mpc}$ box.

Key Contributions

• Hamlet-PM Framework: A novel extension of the Hamlet method that enables direct constraint of ICs from $z=0$ velocity data using a differentiable, non-linear gravity solver, moving beyond the linear theory limitations of the CLUES pipeline.
• Opt-LUM Metric: The introduction of an "optimized Local Universe Model" (opt-LUM). This algorithm iteratively shifts the target position of observed clusters to the center of mass of matched simulated halos. This captures systematic spatial shifts (positional biases) that the standard LUM penalizes, providing a more robust detection rate and mass estimation.
• Standardized Benchmarking: The application of the Velocity Field Olympics (VFO) and LUM metrics to quantitatively compare the new method against existing velocity-based (CLONES, P23) and redshift-based (CSiBORG, Manticore-Local) constrained simulations.

Results

• Large-Scale Structure: The 100 simulations successfully reproduce the observed large-scale structure within the constrained volume, showing persistent filaments and voids. Non-linear structures appear more pronounced than in linear reconstructions.
• Cluster Detection: Using the opt-LUM metric, the detection rate for 12 prominent nearby clusters (e.g., Virgo, Coma, Shapley) significantly improves over the standard LUM. Most clusters achieve detection rates $>80\%$ at $2\sigma$, with many reaching $3\sigma$ and $4\sigma$ thresholds.
• Mass Estimates: The posterior mass distributions for the matched clusters are generally consistent with literature values (X-ray, SZ, and dynamical probes), particularly when using the opt-LUM counterparts.
• Velocity Field Performance: In the VFO benchmark, Hamlet-PM outperforms previous velocity-based pipelines (CLONES) but is outperformed by redshift-based constrained simulations (Manticore-Local, CSiBORG1).
• Systematic Biases: The study identifies a significant over-abundance of massive halos ($>10^{14} h^{-1}M_\odot$) in the constrained volume compared to $\Lambda$CDM predictions and observational estimates. The particle-mesh solver used for inference appears to underestimate velocity amplitudes in dense regions, forcing the initial conditions to deviate from the cosmological model to fit the data, which results in biased resimulations with too many massive clusters.

Significance and Claims
The paper claims that field-level forward modeling of initial conditions can be successfully applied to actual peculiar velocity observations to produce highly constrained cosmological simulations. In its current form, Hamlet-PM is asserted to overtake the quality of the standard CLUES pipeline (WF/RZA/CRs) currently used in the peculiar-velocity community, primarily due to its ability to model non-linearities and complex observational effects within a self-consistent Bayesian framework.

However, the authors are modest regarding the method's current limitations. They acknowledge that systematic biases, specifically the over-production of massive clusters and the underestimation of velocity contrasts in dense regions by the coarse particle-mesh solver, must be addressed. The paper concludes that while the method represents a leap forward for analyzing peculiar velocity data at the field level, it is still outperformed by redshift-based constrained simulations (like Manticore-Local) due to the higher density and precision of redshift data compared to distance-based peculiar velocities. Future improvements are expected to come from refining the observational model (selection functions, cross-calibration) and improving the density/velocity field modeling to probe deeper into the non-linear regime.