Interpreting observations of the Lyman-alpha forest flux power spectrum requires interpolation between a small number of expensive simulations. We present a Gaussian process emulator modelling the 1D flux power spectrum as a function of the amplitude and slope of the small-scale linear matter power spectrum, and the state of the intergalactic medium at the epoch of interest (2 < z < 4). This parameterisation enables the prediction of the flux power spectrum in extended cosmological models that are not explicitly included in the training set, eliminating the need to construct bespoke emulators for a number of extensions to Lambda CDM. Our emulator is appropriate for cosmologies in which the linear matter power spectrum is described to percent level accuracy by just an amplitude and slope across the epoch of interest, and in the regime probed by eBOSS/DESI data. We demonstrate this for massive neutrino cosmologies, where the emulator is able to predict the flux power spectrum in a Sigma m(nu) = 0.3 eV neutrino cosmology to sub-percent accuracy, without including massive neutrinos in the training simulations. Further parameters would be required to describe models with sharp features in the linear power, such as warm or light axion dark matter. This work will facilitate the combination of upcoming DESI data with observations of the cosmic microwave background, to obtain constraints on neutrino mass and other extensions to Lambda CDM cosmology.

We present a general framework for obtaining robust bounds on the nature of dark matter using cosmological N-body simulations and Lyman-alpha forest data. We construct an emulator of hydrodynamical simulations, which is a flexible, accurate and computationally efficient model for predicting the response of the Lyman-alpha forest flux power spectrum to different dark matter models, the state of the intergalactic medium (IGM) and the primordial power spectrum. The emulator combines a flexible parametrization for the small-scale suppression in the matter power spectrum arising in noncold dark matter models, with an improved IGM model. We then demonstrate how to optimize the emulator for the case of ultralight axion dark matter, presenting tests of convergence. We also carry out cross-validation tests of the accuracy of flux power spectrum prediction. This framework can be optimized for the analysis of many other dark matter candidates, e.g., warm or interacting dark matter. Our work demonstrates that a combination of an optimized emulator and cosmological effective theories, where many models are described by a single set of equations, is a powerful approach for robust and computationally efficient inference from the cosmic large-scale structure.

We present a new bound on the ultralight axion (ULA) dark matter mass m(a), using the Lyman-alpha forest to look for suppressed cosmic structure growth: a 95% lower limit m(a) > 2 x 10(-20) eV. This strongly disfavors (>99.7% credibility) the canonical ULA with 10(-22) eV < m(a) < 10(-21) eV, motivated by the string axiverse and solutions to possible tensions in the cold dark matter model. We strengthen previous equivalent bounds by about an order of magnitude. We demonstrate the robustness of our results using an optimized emulator of improved hydrodynamical simulations.

Using a suite of hydrodynamical simulations with cold dark matter, baryons, and neutrinos, we present a detailed study of the effect of massive neutrinos on the 1-D and 3-D flux power spectra of the Lyman-alpha (Ly alpha) forest. The presence of massive neutrinos in cosmology induces a scale- and time-dependent suppression of structure formation that is strongest on small scales. Measuring this suppression is a key method for inferring neutrino masses from cosmological data, and is one of the main goals of ongoing and future surveys like eBOSS, DES, LSST, Euclid or DESI. The clustering in the Ly alpha forest traces the quasilinear power at late times and on small scales. In combination with observations of the cosmic microwave background, the forest therefore provides some of the tightest constraints on the sum of the neutrino masses. However there is a well-known degeneracy between Sigma m(v) and the amplitude of perturbations in the linear matter power spectrum. We study the corresponding degeneracy in the 1-D flux power spectrum of the Ly alpha forest, and for the first time also study this degeneracy in the 3-D flux power spectrum. We show that the non-linear effects of massive neutrinos on the Lycy forest, beyond the effect of linear power amplitude suppression, are negligible, and this degeneracy persists in the Ly alpha forest observables to a high precision. We discuss the implications of this degeneracy for choosing parametrisations of the Ly alpha forest for cosmological analysis.

We present methods for interpolating between the 1-D flux power spectrum of the Lyman-alpha forest, as output by cosmological hydrodynamic simulations. Interpolation is necessary for cosmological parameter estimation due to the limited number of simulations possible. We construct an emulator for the Lyman-alpha forest flux power spectrum from 21 small simulations using Latin hypercube sampling and Gaussian process interpolation. We show that this emulator has a typical accuracy of 1 : 5% and a worst-case accuracy of 4%, which compares well to the current statistical error of 3-5% at z < 3 from BOSS DR9. We compare to the previous state of the art, quadratic polynomial interpolation. The Latin hypercube samples the entire volume of parameter space, while quadratic polynomial emulation samples only lower-dimensional subspaces. The Gaussian process provides an estimate of the emulation error and we show using test simulations that this estimate is reasonable. We construct a likelihood function and use it to show that the posterior constraints generated using the emulator are unbiased. We show that our Gaussian process emulator has lower emulation error than quadratic polynomial interpolation and thus produces tighter posterior confidence intervals, which will be essential for future Lyman-alpha surveys such as DESI.

The Lyman-alpha forest provides strong constraints on both cosmological parameters and intergalactic medium astrophysics, which are forecast to improve further with the next generation of surveys including eBOSS and DESI. As is generic in cosmological inference, extracting this information requires a likelihood to be computed throughout a high-dimensional parameter space. Evaluating the likelihood requires a robust and accurate mapping between the parameters and observables, in this case the 1D flux power spectrum. Cosmological simulations enable such a mapping, but due to computational time constraints can only be evaluated at a handful of sample points; emulators are designed to interpolate between these. The problem then reduces to placing the sample points such that an accurate mapping is obtained while minimising the number of expensive simulations required. To address this, we introduce an emulation procedure that employs Bayesian optimisation of the training set for a Gaussian process interpolation scheme. Starting with a Latin hypercube sampling (other schemes with good space-filling properties can be used), we iteratively augment the training set with extra simulations at new parameter positions which balance the need to reduce interpolation error while focussing on regions of high likelihood. We show that smaller emulator error from the Bayesian optimisation propagates to smaller widths on the posterior distribution. Even with fewer simulations than a Latin hypercube, Bayesian optimisation shrinks the 95% credible volume by 90% and, e.g., the 1 sigma error on the amplitude of small-scale primordial fluctuations by 38%. This is the first demonstration of Bayesian optimisation applied to large-scale structure emulation, and we anticipate the technique will generalise to many other probes such as galaxy clustering, weak lensing and 21cm.

We test a method to reduce unwanted sample variance when predicting Ly alpha forest power spectra from cosmological hydrodynamical simulations. Sample variance arises due to sparse sampling of modes on large scales and propagates to small scales through nonlinear gravitational evolution. To tackle this, we generate initial conditions in which the density perturbation amplitudes are fixed to the ensemble average power spectrum-and are generated in pairs with exactly opposite phases. We run 50 such simulations (25 pairs) and compare their performance against 50 standard simulations by measuring the Ly alpha 1D and 3D power spectra at redshifts z = 2, 3, and 4. Both ensembles use periodic boxes of 40 h(-1)Mpc containing 512(3) particles each of dark matter and gas. As a typical example of improvement, for wavenumbers k = 0.25 hMpc(-1) at z = 3, we find estimates of the 1D and 3D power spectra converge 34 and 12 times faster in a paired-fixed ensemble compared with a standard ensemble. We conclude that, by reducing the computational time required to achieve fixed accuracy on predicted power spectra, the method frees up resources for exploration of varying thermal and cosmological parameters-ultimately allowing the improved precision and accuracy of statistical inference.

Correlations measured in three dimensions in the Lyman-alpha forest are contaminated by the presence of the damping wings of high column density (HCD) absorbing systems of neutral hydrogen (H I; having column densities N(H I) > 1.6 x 10(17) atoms cm(-2)), which extend significantly beyond the redshift-space location of the absorber. We measure this effect as a function of the column density of the HCD absorbers and redshift by measuring threedimensional (3D) flux power spectra in cosmological hydrodynamical simulations from the Illustris project. Survey pipelines exclude regions containing the largest damping wings. We find that, even after this procedure, there is a scale-dependent correction to the 3D Lyman alpha forest flux power spectrum from residual contamination. We model this residual using a simple physical model of the HCD absorbers as linearly biased tracers of the matter density distribution, convolved with their Voigt profiles and integrated over the column density distribution function. We recommend the use of this model over existing models used in data analysis, which approximate the damping wings as top-hats and so miss shape information in the extended wings. The simple 'linear Voigt model' is statistically consistent with our simulation results for a mock residual contamination up to small scales (vertical bar k vertical bar < 1 h Mpc(-1)). It does not account for the effect of the highest column density absorbers on the smallest scales (e.g. vertical bar k vertical bar > 0.4 h Mpc(-1) for small damped Lyman-alpha absorbers; HCD absorbers with N(H-I) similar to 10(21) atoms cm(-2)). However, these systems are in any case preferentially removed from survey data. Our model is appropriate for an accurate analysis of the baryon acoustic oscillations feature. It is additionally essential for reconstructing the full shape of the 3D flux power spectrum.