Second-order Lagrangian Perturbation Theory Research Articles

Large-scale density and velocity field reconstructions with neural networks

ABSTRACT We assess a neural network (NN) method for reconstructing 3D cosmological density and velocity fields (target) from discrete and incomplete galaxy distributions (input). We employ second-order Lagrangian Perturbation Theory to generate a large ensemble of mock data to train an auto-encoder (AE) architecture with a Mean Squared Error (MSE) loss function. The AE successfully captures non-linear features arising from gravitational dynamics and the discreteness of the galaxy distribution. It preserves the positivity of the reconstructed density field and exhibits a weaker suppression of the power on small scales than the traditional linear Wiener filter (WF), which we use as a benchmark. In the density reconstruction, the reduction of the AE MSE relative to the WF is $\sim 15~{{\ \rm per\ cent}}$ , whereas for the velocity reconstruction a relative reduction of up to a factor of two can be achieved. The AE is advantageous to the WF at recovering the distribution of the target fields, especially at the tails. In fact, trained with an MSE loss, any NN estimate approaches the unbiased mean of the underlying target given the input. This implies a slope of unity in the linear regression of the true on the NN-reconstructed field. Only for the special case of Gaussian fields, the NN and WF estimates are equivalent. Nonetheless, we also recover a linear regression slope of unity for the WF with non-Gaussian fields.

Numerical implementation of the Cubic Galileon model in pinocchio

ABSTRACT We present a perturbative treatment of non-linear galaxy clustering in the context of the cubic Galileon modified gravity model, in terms of second-order Lagrangian Perturbation theory and an extension of ellipsoidal collapse that includes Vainshtein screening. We numerically implement such prescriptions in the approximate pinocchio code, and use it to generate realizations of the matter density field and halo catalogues with different prescriptions for ellipsoidal collapse. We investigate the impact of three different approximations in the computation of collapse times on the halo mass function, halo bias, and matter power spectrum. In the halo mass function, both the modified gravity effect and the screening effect are significant in the high-mass end, similar to what is found for other MG models. We perform a comparison with N-body simulations to assess the validity of our approach, and show that we can reproduce the same trend observed in simulations for all quantities considered. With a simple modification to the grouping algorithm of pinocchio to take into account the gravity modification, and without the need to re-calibrate the algorithm, we show that we can reproduce the linear halo bias and the mildly non-linear matter power spectrum of simulations with good accuracy, especially for the implementation with Vainshtein screening. We stress that, while approximate, our method is orders of magnitude faster than a full N-body simulation, making it an optimal tool for the quick generation of large sets of halo catalogues for cosmological observables.

Simulating cosmic structure formation with the gadget-4 code

ABSTRACT Numerical methods have become a powerful tool for research in astrophysics, but their utility depends critically on the availability of suitable simulation codes. This calls for continuous efforts in code development, which is necessitated also by the rapidly evolving technology underlying today’s computing hardware. Here, we discuss recent methodological progress in the gadget code, which has been widely applied in cosmic structure formation over the past two decades. The new version offers improvements in force accuracy, in time-stepping, in adaptivity to a large dynamic range in time-scales, in computational efficiency, and in parallel scalability through a special MPI/shared-memory parallelization and communication strategy, and a more-sophisticated domain decomposition algorithm. A manifestly momentum conserving fast multipole method (FMM) can be employed as an alternative to the one-sided TreePM gravity solver introduced in earlier versions. Two different flavours of smoothed particle hydrodynamics, a classic entropy-conserving formulation and a pressure-based approach, are supported for dealing with gaseous flows. The code is able to cope with very large problem sizes, thus allowing accurate predictions for cosmic structure formation in support of future precision tests of cosmology, and at the same time is well adapted to high dynamic range zoom-calculations with extreme variability of the particle number density in the simulated volume. The gadget-4 code is publicly released to the community and contains infrastructure for on-the-fly group and substructure finding and tracking, as well as merger tree building, a simple model for radiative cooling and star formation, a high dynamic range power spectrum estimator, and an initial condition generator based on second-order Lagrangian perturbation theory.

Transients from initial conditions based on Lagrangian perturbation theory in N-body simulations III: The case of Gadget-2 code

In modern cosmology, the precision of the theoretical prediction is increasingly required. In cosmological [Formula: see text]-body simulations, the effect of higher-order Lagrangian perturbation on the initial conditions appears in terms of statistical quantities of matter density field. We have considered the effect of third-order Lagrangian perturbation theory (3LPT) on the initial conditions, which can be applied to Gadget-2 code. Then, as statistical quantities, non-Gaussianity of matter density field has been compared between cases of different order perturbations for the initial conditions. Then, we demonstrate the validity of the initial conditions with second-order Lagrangian perturbation theory (2LPT).

Fast numerical method to generate halo catalogues in modified gravity (part I): second-order Lagrangian perturbation theory

ABSTRACT We present a new numerical method to determine second-order Lagrangian displacement fields in presence of modified gravity (MG). We start from the extension of Lagrangian perturbation theory (LPT) to a class of MG models, which can be described by a parametrized Poisson equation. We exploit Fast Fourier transforms to compute the full source term of the differential equation for the second-order Lagrangian displacement field. We compare its mean to the source term computed for specific configurations, for which a k-dependent solution can be found numerically. We choose the configuration that best matches the full source term, thus obtaining an approximate factorization of the second-order displacement field as the space term valid for Λ Cold Dark Matter (ΛCDM) times a k-dependent, second-order growth factor. Such approximation is used to compute second-order displacements for particles. The method is tested against N-body simulations run with standard and f(R) gravity: we rely on the results of a friends-of-friends code run on the N-body snapshots to assign particles to haloes, then compute the halo power spectrum. We find very consistent results for the two gravity theories: second-order LPT (2LPT) allows to recover the N-body halo power spectrum within ∼10 per cent precision to k ∼ 0.2–0.4 h Mpc−1, as well as halo positions. We show that the performance of 2LPT with MG is the same (within 1 per cent) as the one obtained for standard ΛCDM case. This formulation of 2LPT can quickly generate dark matter distributions with f(R) gravity, and can easily be extended to other MG theories.

Semiclassical path to cosmic large-scale structure

We chart a path toward solving for the nonlinear gravitational dynamics of cold dark matter by relying on a semiclassical description using the propagator. The evolution of the propagator is given by a Schr\"odinger equation, where the small parameter $\hbar$ acts as a softening scale that regulates singularities at shell-crossing. The leading-order propagator, called free propagator, is the semiclassical equivalent of the Zel'dovich approximation (ZA), that describes inertial particle motion along straight trajectories. At next-to-leading order, we solve for the propagator perturbatively and obtain, in the classical limit the displacement field from second-order Lagrangian perturbation theory (LPT). The associated velocity naturally includes an additional term that would be considered as third order in LPT. We show that this term is actually needed to preserve the underlying Hamiltonian structure, and ignoring it could lead to the spurious excitation of vorticity in certain implementations of second-order LPT. We show that for sufficiently small $\hbar$ the corresponding propagator solutions closely resemble LPT, with the additions that spurious vorticity is avoided and the dynamics at shell-crossing is regularised. Our analytical results possess a symplectic structure that allows us to advance numerical schemes for the large-scale structure. For times shortly after shell-crossing, we explore the generation of vorticity, which in our method does not involve any explicit multi-stream averaging, but instead arises naturally as a conserved topological charge.

Grid-based calculation for perturbation theory of large-scale structure

Perturbation theory (PT) calculation of large-scale structure has been used to interpret the observed non-linear statistics of large-scale structure at the quasi-linear regime. In particular, the so-called standard perturbation theory (SPT) provides a basis for the analytical computation of the higher-order quantities of large-scale structure. Here, we present a novel, grid-based algorithm for the SPT calculation, hence named GridSPT, to generate the higher-order density and velocity fields from a given linear power spectrum. Taking advantage of the Fast Fourier Transform, the GridSPT quickly generates the nonlinear density fields at each order, from which we calculate the statistical quantities such as non-linear power spectrum and bispectrum. Comparing the density fields (to fifth order) from GridSPT with those from the full N-body simulations in the configuration space, we find that GridSPT accurately reproduces the N-body result on large scales. The agreement worsens with smaller smoothing radius, particularly for the under-dense regions where we find that 2LPT (second-order Lagrangian perturbation theory) algorithm performs well.

Toward accurate modelling of the non-linear matter bispectrum: standard perturbation theory and transients from initial conditions

Accurate modeling of nonlinearities in the galaxy bispectrum, the Fourier transform of the galaxy three-point correlation function, is essential to fully exploit it as a cosmological probe. In this paper, we present numerical and theoretical challenges in modeling the nonlinear bispectrum. First, we test the robustness of the matter bispectrum measured from N-body simulations using different initial conditions generators. We run a suite of N-body simulations using the Zel'dovich approximation and second-order Lagrangian perturbation theory (2LPT) at different starting redshifts, and find that transients from initial decaying modes systematically reduce the nonlinearities in the matter bispectrum. To achieve 1% accuracy in the matter bispectrum for $z\le3$ on scales $k<1$ $h$/Mpc, 2LPT initial conditions generator with initial redshift of $z\gtrsim100$ is required. We then compare various analytical formulas and empirical fitting functions for modeling the nonlinear matter bispectrum, and discuss the regimes for which each is valid. We find that the next-to-leading order (one-loop) correction from standard perturbation theory matches with N-body results on quasi-linear scales for $z\ge1$. We find that the fitting formula in Gil-Mar\'{\i}n et al. (2012) accurately predicts the matter bispectrum for $z\le1$ on a wide range of scales, but at higher redshifts, the fitting formula given in Scoccimarro & Couchman (2001) gives the best agreement with measurements from N-body simulations.

Truthing the stretch: non-perturbative cosmological realizations with multiscale spherical collapse

Abstract Here we present a simple, parameter-free, non-perturbative algorithm that gives low-redshift cosmological particle realizations accurate to few-Megaparsec scales, called muscle (MUltiscale Spherical-ColLapse Evolution). It has virtually the same cost as producing N-body-simulation initial conditions, since it works with the ‘stretch’ parameter ψ, the Lagrangian divergence of the displacement field. It promises to be useful in quickly producing mock catalogues, and to simplify computationally intensive reconstructions of galaxy surveys. muscle applies a spherical-collapse prescription on multiple Gaussian-smoothed scales. It achieves higher accuracy than perturbative schemes (Zel'dovich and second-order Lagrangian perturbation theory – 2LPT), and, by including the void-in-cloud process (voids in large-scale collapsing regions), solves problems with a single-scale spherical-collapse scheme. Slight further improvement is possible by mixing in the 2LPT estimate on large scales. Additionally, we show the behaviour of ψ for different morphologies (voids, walls, filaments, and haloes). A python code to produce these realizations is available at http://skysrv.pha.jhu.edu/∼neyrinck/muscle.html.

EUNHA: A NEW COSMOLOGICAL HYDRODYNAMIC SIMULATION CODE

We develop a parallel cosmological hydrodynamic simulation code designed for the study of formation and evolution of cosmological structures. The gravitational force is calculated using the TreePM method and the hydrodynamics is implemented based on the smoothed particle hydrodynamics. The initial displacement and velocity of simulation particles are calculated according to second-order Lagrangian perturbation theory using the power spectra of dark matter and baryonic matter. The initial background temperature is given by Recfast and the temperature fluctuations at the initial particle position are assigned according to the adiabatic model. We use a time-limiter scheme over the individual time steps to capture shock-fronts and to ease the time-step tension between the shock and preshock particles. We also include the astrophysical gas processes of radiative heating/cooling, star formation, metal enrichment, and supernova feedback. We test the code in several standard cases such as one-dimensional Riemann problems, Kelvin-Helmholtz, and Sedov blast wave instability. Star formation on the galactic disk is investigated to check whether the Schmidt-Kennicutt relation is properly recovered. We also study global star formation history at different simulation resolutions and compare them with observations.

The initial conditions of the Universe from constrained simulations

Abstract I present a new approach to recover the primordial density fluctuations and the cosmic web structure underlying a galaxy distribution. The method is based on sampling Gaussian fields which are compatible with a galaxy distribution and a structure formation model. This is achieved by splitting the inversion problem into two Gibbs-sampling steps: the first being a Gaussianization step transforming a distribution of point sources at Lagrangian positions – which are not a priori given– into a linear alias-free Gaussian field. This step is based on Hamiltonian sampling with a Gaussian–Poisson model. The second step consists on a likelihood comparison in which the set of matter tracers at the initial conditions is constrained on the galaxy distribution and the assumed structure formation model. For computational reasons second-order Lagrangian perturbation theory is used. However, the presented approach is flexible to adopt any structure formation model. A semi-analytic halo-model-based galaxy mock catalogue is taken to demonstrate that the recovered initial conditions are closely unbiased with respect to the actual ones from the corresponding N-body simulation down to scales of a ∼5 Mpc h−1. The cross-correlation between them shows a substantial gain of information, being at k ∼ 0.3 h Mpc−1 more than doubled. In addition the initial conditions are extremely well Gaussian distributed and the power spectra follow the shape of the linear power spectrum being very close to the actual one from the simulation down to scales of k ∼ 1 h Mpc−1.

Cosmic structure and dynamics of the local Universe

Abstract We present a cosmography analysis of the local Universe based on the recently released Two-Micron All-Sky Redshift Survey catalogue. Our method is based on a Bayesian Networks Machine Learning algorithm (the Kigen-code) which self-consistently samples the initial density fluctuations compatible with the observed galaxy distribution and a structure formation model given by second-order Lagrangian perturbation theory (2LPT). From the initial conditions we obtain an ensemble of reconstructed density and peculiar velocity fields which characterize the local cosmic structure with high accuracy unveiling non-linear structures like filaments and voids in detail. Coherent redshift-space distortions are consistently corrected within 2LPT. From the ensemble of cross-correlations between the reconstructions and the galaxy field and the variance of the recovered density fields, we find that our method is extremely accurate up to k∼1 h Mpc−1 and still yields reliable results down to scales of about 3–4 h−1 Mpc. The motion of the Local Group we obtain within ∼80 h−1 Mpc (νLG = 522 ± 86 km s−1, lLG = 291° ± 16°, bLG = 34° ± 8°) is in good agreement with measurements derived from the cosmic microwave background and from direct observations of peculiar motions and is consistent with the predictions of ΛCDM.

Estimating cosmic velocity fields from density fields and tidal tensors

In this work we investigate the non-linear and non-local relation between cosmological density and peculiar velocity fields. Our goal is to provide an algorithm for the reconstruction of the non-linear velocity field from the fully non-linear density. We find that including the gravitational tidal field tensor using second-order Lagrangian perturbation theory based upon an estimate of the linear component of the non-linear density field significantly improves the estimate of the cosmic flow in comparison to linear theory not only in the low density, but also and more dramatically in the high-density regions. In particular we test two estimates of the linear component: the lognormal model and the iterative Lagrangian linearization. The present approach relies on a rigorous higher order Lagrangian perturbation theory analysis which incorporates a non-local relation. It does not require additional fitting from simulations being in this sense parameter free, it is independent of statistical–geometrical optimization and it is straightforward and efficient to compute. The method is demonstrated to yield an unbiased estimator of the velocity field on scales ≳5 h−1 Mpc with closely Gaussian distributed errors. Moreover, the statistics of the divergence of the peculiar velocity field is extremely well recovered showing a good agreement with the true one from N-body simulations. The typical errors of about 10 km s−1 (1σ confidence intervals) are reduced by more than 80 per cent with respect to linear theory in the scale range between 5 and 10 h−1 Mpc in high-density regions (δ > 2). We also find that iterative Lagrangian linearization is significantly superior in the low-density regime with respect to the lognormal model.

Second-order Lagrangian perturbation theory initial conditions for resimulations

We describe and test a new method for creating initial conditions for cosmological N-body dark matter simulations based on second-order Lagrangian perturbation theory (2lpt). The method can be applied to multi-mass particle distributions making it suitable for creating resimulation, or `zoom' initial conditions. By testing against an analytic solution we show that the method works well for a spherically symmetric perturbation with radial features ranging over more than three orders of magnitude in linear scale and eleven orders of magnitude in particle mass. We apply the method and repeat resimulations of the rapid formation of a high mass halo at redshift fifty and the formation of a Milky-Way mass dark matter halo at redshift zero. In both cases we find that many properties of the final halo show a much smaller sensitivity to the start redshift with the 2lpt initial conditions, than simulations started from Zel'dovich initial conditions. For spherical overdense regions structure formation is erroneously delayed in simulations starting from Zel'dovich initial conditions, and we demonstrate for the case of the formation the high redshift halo that this delay can be accounted for using the spherical collapse model. In addition to being more accurate, 2lpt initial conditions allow simulations to start later, saving computer time.

Transients from initial conditions in cosmological simulations

We study the impact of setting initial conditions in numerical simulations using the standard procedure based on the Zel'dovich approximation (ZA). As it is well known from the perturbation theory, ZA initial conditions have incorrect second- and higher-order growth and therefore excite long-lived transients in the evolution of the statistical properties of density and velocity fields. We also study the improvement brought by using more accurate initial conditions based on second-order Lagrangian perturbation theory (2LPT). We show that 2LPT initial conditions reduce transients significantly and thus are much more appropriate for numerical simulations devoted to precision cosmology. Using controlled numerical experiments with ZA and 2LPT initial conditions, we show that simulations started at redshift z i = 49 using the ZA underestimate the power spectrum in the non-linear regime by about 2,4 and 8 per cent at z = 0, 1, and 3, respectively, whereas the mass function of dark matter haloes is underestimated by 5 per cent at m = 10 15 M ⊙ h -1 (z = 0) and 10 per cent at m = 2 × 10 14 M ⊙ h - (z = 1). The clustering of haloes is also affected to the few per cent level at z = 0. These systematics effects are typically larger than statistical uncertainties in recent mass function and power spectrum fitting formulae extracted from numerical simulations. At large scales, the measured transients in higher-order correlations can be understood from first principle calculations based on perturbation theory.

Higher Order Contributions to the 21 cm Power Spectrum

We consider the contribution of third- and fourth-order terms to the 21 cm power spectrum during the epoch of reionization, which arise because the 21 cm brightness temperature involves a product of the hydrogenic neutral fraction and the gas density. The third-order terms vanish for Gaussian random fields and have been previously neglected or ignored. We measure these terms from radiative transfer simulations and estimate them using cosmological perturbation theory. In our simulated models, the higher order terms are significant, and neglecting them leads to inaccurate 21 cm power spectrum estimates. On small scales the higher order terms are produced by gravitational mode coupling. Small-scale structure grows more readily in large-scale overdense regions, but the same regions tend to be ionized and hence do not contribute to the 21 cm signal. This modifies an earlier intuition that the 21 cm power spectrum simply traces the density power spectrum on scales smaller than that of a typical bubble and implies that small-scale measurements contain more information about the ionizing sources than previously believed. On large scales, higher order moments are not directly related to gravity. They are nonzero because overdense regions tend to ionize first and are important in magnitude at late times, due to large fluctuations in the neutral fraction. Finally, we show that second-order Lagrangian perturbation theory approximately reproduces the statistics of the density field from full numerical simulations for all redshifts and scales of interest, including the mode-coupling effects mentioned above. It can, therefore, be used in conjunction with semianalytic models to explore the broad regions of parameter space relevant for future 21 cm surveys.

The Bispectrum: From Theory to Observations

The bispectrum is the lowest order statistic sensitive to the shape of structures generated by gravitational instability and is a potentially powerful probe of galaxy biasing and the Gaussianity of primordial fluctuations. Although the evolution of the bispectrum is well understood theoretically from nonlinear perturbation theory and numerical simulations, applications to galaxy surveys require a number of issues to be addressed. In this paper we consider the effect on the bispectrum of stochastic nonlinear biasing, radial redshift distortions, non-Gaussian initial conditions, survey geometry, and sampling. We find that: (1) Bias stochasticity does not affect the use of the bispectrum to recover the mean biasing relation between galaxies and mass, at least for models in which the scatter is uncorrelated at large scales. (2) Radial redshift distortions do not significantly change the monopole power spectrum and bispectrum compared to their plane-parallel values. (3) Survey geometry leads to finite-volume effects, which must be taken into account in current surveys before comparison with theoretical predictions can be made. (4) Sparse sampling and survey geometry correlate different triangles, leading to a breakdown of the Gaussian likelihood approximation. We develop a likelihood analysis using bispectrum eigenmodes, calculated by Monte Carlo realizations of mock surveys generated with second-order Lagrangian perturbation theory and checked against N-body simulations. In a companion paper, we apply these results to the analysis of the bispectrum of IRAS galaxies.

Non-linear evolution of the tidal angular momentum of protostructures — II. Non-Gaussian initial conditions

The formalism that describes the non-linear growth of the angular momentum L of protostructures from tidal torques in a Friedmann Universe, as developed in a previous paper, is extended to include non-Gaussian initial conditions. We restrict our analysis here to a particular class of non-Gaussian primordial distributions, namely multiplicative models. In such models, strongly correlated phases are produced by obtaining the gravitational potential via a nonlinear local transformation of an underlying Gaussian random field. The dynamical evolution of the system is followed by describing the trajectories of fluid particles using second-order Lagrangian perturbation theory. In the Einstein-de Sitter universe, the lowest-order perturbative correction to the variance of the linear angular momentum of collapsing structures grows as t^8/3 for generic non-Gaussian statistics, which contrasts with the t^10/3 growth rate characteristic of Gaussian statistics. This is a consequence of the fact that the lowest-order perturbative spin contribution in the non-Gaussian case arises from the third moment of the gravitational potential, which is identically zero for a Gaussian field. Evaluating these corrections at the maximum expansion time of the collapsing structure, we find that these non-Gaussian and non-linear terms can be as high as the linear estimate, without the degree of non-Gaussianity as quantified by skewness and kurtosis of the density field being unacceptably large. The results suggest that higher-order terms in the perturbative expansion may contribute significantly to galactic spin which contrasts with the straightforward Gaussian case.

Second-order solution for determining density and velocity fields of galaxies

In this Letter, I use second-order Lagrangian perturbation theory to calculate an analytical expression relating density to velocity in a gravitating system. This solution can be used to compare peculiar velocity field measurements with observations of large-scale structure. The predictions of both linear theory and second-order theory are compared with the results of N-body simulations. While linear theory systematically overestimates the velocity flows in high-density regions, the second-order corrections calculated herein remove this systematic error.