Einstein’s Universe: Cosmological structure formation in numerical relativity

We perform large-scale cosmological simulations that solve Einstein’s equations directly via numerical relativity. Starting with initial conditions sampled from the cosmic microwave background, we track the emergence of a cosmic web without the need for a background cosmology. We measure the backreaction of large-scale structure on the evolution of averaged quantities in a matter-dominated universe. Although our results are preliminary, we find the global backreaction energy density is of order ${10}^{-8}$ compared to the energy density of matter in our simulations, and is thus unlikely to explain accelerating expansion under our assumptions. Sampling scales above the homogeneity scale of the Universe ($100–180h^{-1}\mathrm{Mpc}$), in our chosen gauge, we find 2–3% variations in local spatial curvature.

Figure 1

Matter power spectrum of our initial conditions. The grey dashed curve shows the power spectrum produced with camb. We show the power as a function of wave number $|\mathbf{k}|=\sqrt{k_x^2+k_y^2+k_z^2}$. The magenta curve shows the section of the power spectrum we sample when using a domain size of $L=1\mathrm{Gpc}$ with resolution $256^3$.

Figure 2

Initial conditions drawn from the cosmic microwave background power spectrum. Here we show initial conditions for the density (top row), metric (middle row), and velocity (bottom row) perturbations for three different physical domain sizes. Left to right shows domain sizes $L=1\mathrm{Gpc}$, 500 Mpc, and 100 Mpc. We show a two-dimensional slice through the midplane of each domain. All initial conditions shown here are at $256^3$ resolution, and all quantities are shown in code units, normalized by the speed of light for the metric and velocity perturbations. The magnitude of the velocity is $|v|=\sqrt{v_x^2+v_y^2+v_z^2}$.

Figure 3

Evolution of a fully general-relativistic cosmic web. Here we show a $256^3$ simulation, in an $L=1\mathrm{Gpc}$ domain. This simulation has evolved from the cosmic microwave background ($z=1100$; top left) until today ($z=0$; bottom right). Each panel shows a two-dimensional slice of the density perturbation in the midplane of the domain. We can see the familiar web structure of modern cosmological $N$-body simulations using Newtonian gravity; however this cosmic web contains all of the corresponding general-relativistic information. The standard deviations of the fractional density perturbation $\delta$ for each panel (progressing in time) are ${\sigma }_{\delta }=0.0026$, 0.15, 0.6, 1.11, 1.89, and 3.92, respectively.

Figure 4

Scale dependence of the cosmic web. The three separate simulations were computed at a resolution of $128^3$ (left to right) with domain sizes $L=1\mathrm{Gpc}$, 500 Mpc, and 100 Mpc, respectively. All snapshots show a two-dimensional density slice in the midplane of the simulation domain at redshift $z=0$.

Figure 5

General-relativistic attributes of an inhomogeneous, anisotropic universe. Panels (left to right) show the matter expansion rate $\theta$, the spatial Ricci curvature $R$, and the shear ${\sigma }^2$, respectively, each relative to the global Hubble expansion ${\mathcal{H}}_{\mathrm{all}}$. Each panel shows a two-dimensional slice at $z=0$ through the midplane of the $L=1\mathrm{Gpc}$ domain at $256^3$ resolution.

Figure 6

Globally, our expansion coincides with that of FLRW. The blue curve in the top panel shows the effective scale factor $a_{\mathcal{D}}$, calculated over the entire $L=1\mathrm{Gpc}$ domain. The dashed magenta curve shows the equivalent FLRW solution (with ${\mathrm{\Omega }}_m=1$), as a function of redshift. The bottom panel shows the residual error for this $256^3$ resolution calculation.

Figure 7

Growing inhomogeneity in matter, curvature, and backreaction. Here we show the cosmological parameters for spheres with various radii $r_{\mathcal{D}}$, randomly placed within an $L=1\mathrm{Gpc}$ domain at $256^3$ resolution. Black points show mean values over 1000 spheres at each radius, while progressively lighter blue and purple shaded regions show the 68%, 95%, and 99.7% confidence intervals for ${\mathrm{\Omega }}_m$ and ${\mathrm{\Omega }}_R$, respectively. Crosses show the mean contribution from backreaction terms ${\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L$, while dashed, dot-dashed, and dotted lines show the 68%, 95%, and 99.7% confidence intervals, respectively. Left to right panels are redshifts $z=9.9$, 1.1, and 0.0, respectively.

Figure 8

We approach homogeneity when averaging over larger scales. Here we show the right-most panel of Fig. 7 extending to averaging radius $r_{\mathcal{D}}=250\mathrm{Mpc}$. Black points show the mean ${\mathrm{\Omega }}_m$, ${\mathrm{\Omega }}_R$, and ${\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L$ over the 1000 spheres at each radius. Progressively lighter blue and purple shaded regions show the 68%, 95%, and 99.7% confidence intervals for ${\mathrm{\Omega }}_m$ and ${\mathrm{\Omega }}_R$, while dashed, dot-dashed, and dotted lines show these for ${\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L$.

Figure 9

Inhomogeneous expansion as a function of time, showing the effective scale factor $a_{\mathcal{D}}$ calculated in spheres of radius 100 Mpc as a function of global expansion $a_{\mathcal{D},\mathrm{all}}$. We calculate $a_{\mathcal{D}}$ in an $L=500\mathrm{Mpc}$ simulation at $128^3$ resolution. Blue curves show overdense regions with $\delta \ge 0.1$, while purple curves show underdense regions with $\delta \le -0.1$. The black dashed line shows the mean expansion over the whole domain.

Figure 10

Relation between the fractional density perturbation, $\delta$, and the deviation in the Hubble parameter, $\delta {\mathcal{H}}_{\mathcal{D}}/{\mathcal{H}}_{\mathrm{all}}$, for averaging radii $r_{\mathcal{D}}=20$, 40, and 80 Mpc (left to right), respectively. Points in each panel represent individual spheres of 1000 sampled at each radius, and the solid line of the same color is the best-fit linear relation, with the slope indicated in each panel. The black dashed line is the prediction from linear theory.

Figure 11

Relative (right panels) and raw (left panels) constraint violation as a function of effective redshift calculated using Eqs. (c4) and (c2), respectively. We show violations for the simulations with a controlled number of physical modes. The top panels show the Hamiltonian constraint violation, and the bottom panels show the momentum constraint violation. Colors show different resolutions as indicated by the legend.

Figure 12

Second-order convergence of the Hamiltonian and momentum constraints for the set of simulations with a fixed number of physical modes. We show the $L_1$ error calculated using Eq. (c2) for both constraints. The top two panels show the $L_1$ error for the initial conditions, at $z=1100$, and the bottom two panels show the $L_1$ error for $z=0$.

Figure 13

Power spectrum of fractional density fluctuations $\delta$ at $z=0$. Solid colored curves show $P(k)$ for three simulations with $64^3$, $128^3$, $256^3$, each for an $L=1\mathrm{Gpc}$ domain, and the dashed curve shows the linear power spectrum at $z=0$. The pink shaded region represents one-dimensional scales of 80–100 Mpc.

Figure 14

Global cosmological parameters as a function of effective redshift, for simulations with full power spectrum sampling (left panel) and a controlled number of physical modes (right panel). Dotted, dashed, and solid curves show resolutions as indicated in each separate legend. Blue curves show ${\mathrm{\Omega }}_m$, green curves show $|{\mathrm{\Omega }}_R|$, and purple curves show $|{\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L|$.

Figure 15

Cosmological parameters measured at $z=0$ for simulations with a controlled number of physical modes. Colored points show ${\mathrm{\Omega }}_m$ (left), $|{\mathrm{\Omega }}_R|$ (middle), and $|{\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L|$ (right) for resolutions $N=32$, 64, and 128. Dashed curves are the convergence fit for each parameter, detailed in the text. We use these curves for a Richardson extrapolation to calculate the true values of the parameters and hence the errors on our measurements.

Figure 16

Richardson extrapolated errors for cosmological parameters calculated within subdomains, for the $128^3$ resolution controlled simulation, as a function of averaging radius $r_{\mathcal{D}}$. Blue points show the percentage error for ${\mathrm{\Omega }}_m$, green points for ${\mathrm{\Omega }}_R$, and purple points for ${\mathrm{\Omega }}_Q+{\mathrm{\Omega }}_L$.