Zel'dovich Approximation

The Zel'dovich Approximation is a kinematic first-order Lagrangian perturbation theory of cosmological large-scale structure formation formulated by the Soviet astrophysicist Yakov Borisovich Zel'dovich in 1970. Whereas standard Eulerian perturbation theory expands the matter density contrast δ(x, t) = (ρ(x, t) − ρ̄(t)) / ρ̄(t) at fixed spatial coordinates and breaks down as soon as density fluctuations approach unity (|δ| ~ 1), the Zel'dovich approximation tracks the ballistic trajectories of collisionless fluid elements from their initial comoving Lagrangian positions q to their Eulerian comoving positions x(q, t) = q − D(t)∇_qΨ(q). By relating the local nonlinear density to the three eigenvalues (α, β, γ) of the initial deformation tensor Ψ_ij = ∂²Ψ / (∂q_i ∂q_j), Zel'dovich demonstrated that generic three-dimensional gravitational collapse is inherently anisotropic—flattening first along the axis of largest eigenvalue into two-dimensional sheets ("Zel'dovich pancakes"), then intersecting along one-dimensional filaments, and finally draining into virialized cluster nodes to form the cosmic web.

Consider a universe expanding with scale factor a(t), relating proper physical coordinates r to comoving Eulerian coordinates x via r = a(t)x. Each collisionless dark matter or baryonic fluid parcel is labeled by its initial unperturbed comoving coordinate q = x(t → 0) in the nearly homogeneous early universe.

In the Zel'dovich Approximation, the comoving Eulerian position x(q, t) of a fluid element at cosmic time t is separated into a time-dependent growth factor D(t) and a spatial displacement vector s(q) = −∇_qΨ(q) derived from the initial velocity/gravitational potential Ψ(q):

x(q, t) = q + D(t) s(q) = q − D(t) ∇_qΨ(q)

Here D(t) is the standard linear growing mode of gravitational instability (satisfying d²D / dt² + 2H(t)(dD / dt) = 4*πGρ̄*(t)D, so D(t) ∝ a(t) in an Einstein–de Sitter matter-dominated universe), and Ψ(q) is proportional to the primordial peculiar gravitational potential ϕ₀(q) via the Poisson equation:

∇_q²Ψ(q) = δ_L(q, t) / D(t) = δ₀(q)

Taking the time derivative of a(t)x(q, t) yields the peculiar velocity v_pec = a(*d*x / dt) = −a(dD / dt)∇_qΨ(q). In comoving velocity coordinates u ≡ *d*x / dD = −∇_qΨ(q), each fluid element moves at constant velocity *d*u / dD = 0 along a straight ballistic line determined by the initial potential gradient until trajectories cross (shell crossing).

Complete One-Dimensional Solution Before Shell Crossing

A remarkable mathematical property of the Zel'dovich approximation in one spatial dimension (planar symmetry) is that it represents the complete nonlinear solution to the coupled Euler–Poisson equations of collisionless self-gravitating matter all the way up to the moment of shell crossing, because in one dimension the gravitational force exerted on a planar mass sheet is proportional to the column mass enclosed between it and the origin, which remains conserved before sheets cross.

Conservation of mass between the initial uniform Lagrangian volume element d³q (with mean comoving density ρ̄₀) and the deformed Eulerian volume element d³x (with nonlinear density ρ(x, t)) requires:

ρ(x, t) d³x = ρ̄₀ d³q ⇒ 1 + δ(x, t) = ρ(x, t) / ρ̄₀ = | det(∂x_i / ∂q_j) |⁻¹

From x(q, t) = q − D(t)∇_qΨ(q), the Jacobian matrix J_ij = ∂x_i / ∂q_j is symmetric and given by:

J_ij(q, t) = δ_ij − D(t) (∂²Ψ / ∂q_i ∂q_j) ≡ δ_ij − D(t) Ψ_ij(q)

Because the initial strain (deformation) tensor Ψ_ij(q) = ∂²Ψ / (∂q_i ∂q_j) is a real symmetric 3 × 3 matrix, it can be diagonalized locally at every Lagrangian point q along its three orthogonal principal axes with real eigenvalues ordered as α(q) ≥ β(q) ≥ γ(q). In the principal-axis frame, the Jacobian determinant factors into a product over the three principal axes, yielding Zel'dovich's nonlinear density formula:

1 + δ(q, t) = 1 / [ (1 − D(t)α(q)) (1 − D(t)β(q)) (1 − D(t)γ(q)) ]

Recovery of Linear Theory and First-Collapse Singularity

When D(t)|α| ≪ 1, expanding the denominator to first order recovers standard Eulerian linear perturbation theory:

δ_L(q, t) ≈ D(t) (α + β + γ) = D(t) Tr(Ψ_ij) = D(t) ∇_q²Ψ(q)

In a Gaussian random field of primordial curvature perturbations (analyzed statistically by Andrey Doroshkevich in 1970), the probability that all three eigenvalues are simultaneously equal (α = β = γ, corresponding to perfectly spherical collapse as assumed in the Press–Schechter formalism) has measure zero. Instead, generic collapse proceeds sequentially through a three-stage topological hierarchy as D(t) grows:

1. Sheets ("Zel'dovich Pancakes", D(t) = 1 / α): When D(t) reaches 1 / α (while 1 − Dβ > 0 and 1 − Dγ > 0), the first factor (1 − D(t)α) vanishes and the local density formally diverges (ρ → ∞) along a two-dimensional caustic sheet perpendicular to the α̂ eigenvector.
2. Cosmic Filaments (D(t) = 1 / β): If the second eigenvalue β > 0, subsequent collapse along the second principal axis compresses the sheet into a dense one-dimensional filament along the intersection of pancakes.
3. Cluster Nodes (D(t) = 1 / γ): Where all three eigenvalues are positive (α ≥ β ≥ γ > 0), matter drains along filaments into three-dimensionally collapsed galaxy cluster nodes, whereas regions with α < β < γ < 0 expand along all three axes to form cosmic voids.

Once fluid trajectories cross inside a Zel'dovich pancake at D(t) > 1 / α(q), multistreaming occurs. Because the basic Zel'dovich approximation assigns each particle a constant comoving velocity u = −∇_qΨ(q) independent of the evolving local gravitational potential after shell crossing, particles pass straight through the pancake without decelerating, causing the sheet to thicken indefinitely rather than remaining gravitationally bound.

The Burgers Adhesion Model (Gurbatov, Saichev & Shandarin 1989)

To prevent artificial post-caustic dispersal while preserving analytical tractability, Sergei Gurbatov, Alexander Saichev, and Sergei Shandarin introduced the Adhesion Approximation, adding an artificial kinematic viscosity term ν∇_x²u in the vanishing-viscosity limit (ν → 0⁺) to mimic the gravitational sticking of colliding dark matter streams:

∂u / ∂D + (u · ∇_x)u = ν ∇_x²u

This is the three-dimensional potential Burgers equation, which can be solved exactly via the Hopf–Cole transformation u = −2ν∇_x ln U, reducing non-linear cosmic-web skeleton construction to a geometric paraboloid-tangent convex hull construction on the initial potential Ψ(q).

Second-Order Lagrangian Perturbation Theory (2LPT) and BAO Reconstruction

Extending the Lagrangian perturbation expansion to second order (2LPT, Buchert 1992; Bouchet et al. 1995) adds a second-order displacement potential Ψ⁽²⁾(q) sourced by the second principal invariant of Ψ_ij:

x(q, t) = q − D₁(t) ∇_qΨ⁽¹⁾(q) + D₂(t) ∇_qΨ⁽²⁾(q), D₂ ≈ −(3 / 7) Ω_m^−1/143 D₁²

Today, the Zel'dovich approximation and 2LPT serve two foundational roles in precision observational cosmology:

1. Cosmological N-Body Simulation Initial Conditions: Every modern cosmological simulation (such as Millennium, IllustrisTNG, and Euclid Flagship) uses the Zel'dovich approximation or 2LPT at high redshift (z ≈ 49–127) to displace particles from a uniform grid or glass configuration onto growing-mode trajectories without exciting spurious decaying modes.
2. Baryon Acoustic Oscillation (BAO) Density-Field Reconstruction: Proposed by Daniel Eisenstein and collaborators in 2007 and applied across BOSS, eBOSS, and DESI galaxy surveys, BAO reconstruction estimates the large-scale smoothed displacement field s(x) = −∇Ψ(x) from the observed galaxy density field and shifts galaxies backward by −s to undo nonlinear Zel'dovich bulk flows, sharpening the 150-Mpc baryon acoustic oscillations peak and improving distance-scale precision by a factor of ~1.5–2.

Why did Zel'dovich originally develop his pancake theory for hot dark matter and baryons?

When Zel'dovich formulated the approximation in 1970, cosmological models focused on adiabatic baryon perturbations (with Silk damping erasing scales below ~10¹³ M_☉) or massive neutrino hot dark matter (with free-streaming erasing galaxy-scale fluctuations below ~10¹⁴–10¹⁵ M_☉). In such "top-down" scenarios with a sharp small-scale cutoff in the initial power spectrum P(k), supercluster-scale Zel'dovich pancakes are the very first structures to collapse in the universe, subsequently fragmenting into galaxies via radiative shock cooling. In modern cold dark matter (ΛCDM) cosmology, small-scale halos collapse hierarchically ("bottom-up"), but the Zel'dovich approximation applied to the smoothed large-scale potential still governs the megaparsec-scale cosmic web of sheets, filaments, and voids.

How does the Zel'dovich approximation relate to the Press–Schechter halo mass function?

The classic 1974 Press–Schechter formalism assumes isotropic spherical collapse, in which a region collapses when its linearly extrapolated trace overdensity δ_L = D(t)(α + β + γ) exceeds the spherical threshold δ_c ≈ 1.686. By incorporating the anisotropic tidal shear eigenvalues (α, β, γ) of the Zel'dovich deformation tensor into an ellipsoidal collapse barrier (where shear delays small-halo collapse and accelerates high-mass collapse), Ravi Sheth, H. J. Mo, and Giuseppe Tormen (2001) derived the Sheth–Tormen halo mass function, resolving the systematic discrepancies of spherical Press–Schechter against N-body simulations.

What happens to baryonic gas when a Zel'dovich pancake collapses?

Unlike collisionless dark matter particles—which pass through the midplane caustic to form a multi-stream phase-space sheet—collisional baryonic gas cannot cross itself. Instead, as infalling gas streams collide at the midplane of a Zel'dovich pancake, they form a pair of outward-propagating accretion shocks that convert the supersonically infalling kinetic energy into thermal X-ray-emitting heat (T ~ 10⁶–10⁷ K in the warm–hot intergalactic medium, or WHIM), trapping a dense, shock-heated gas layer between the two shock fronts.

1. Zel'dovich, Y. B. (1970). Gravitational instability: An approximate theory for large density perturbations. Astronomy and Astrophysics, 5, 84–89. https://ui.adsabs.harvard.edu/abs/1970A%26A.....5...84Z
2. Doroshkevich, A. G. (1970). Spatial structure of perturbations and origin of galactic rotation in fluctuation theory. Astrophysics (Astrofizika), 6(4), 320–330. https://doi.org/10.1007/BF01001625
3. Shandarin, S. F., & Zel'dovich, Y. B. (1989). The large-scale structure of the universe: Turbulence, intermittency, structures in a self-gravitating medium. Reviews of Modern Physics, 61(2), 185–220. https://doi.org/10.1103/RevModPhys.61.185
4. Bouchet, F. R., Colombi, S., Hivon, E., & Juszkiewicz, R. (1995). Perturbative Lagrangian approach to gravitational instability. Astronomy and Astrophysics, 296, 575–608. https://ui.adsabs.harvard.edu/abs/1995A%26A...296..575B
5. Eisenstein, D. J., Seo, H.-J., Sirko, E., & Spergel, D. N. (2007). Improving cosmological distance measurements by reconstruction of the baryon acoustic peak. The Astrophysical Journal, 664(2), 675–679. https://doi.org/10.1086/518712