Zel'dovich Approximation

A universe that collapses into sheets

Roughly 380,000 years after the Big Bang, the matter in the universe was almost perfectly smooth — density variations of only about one part in 100,000. Today that same matter is concentrated into a frothy lattice of walls, filaments, and dense knots separated by enormous near-empty voids: the cosmic web. The question of how you get from one to the other is the central problem of structure formation, and the first theory to predict the shape of the answer — not just that clumps grow, but that they grow into sheets — was Yakov Zel'dovich's 1970 approximation.

The Zel'dovich approximation is deceptively simple. Take every mass element, label it by its initial position q (its "Lagrangian coordinate"), and let it coast in a straight line along a fixed direction set by the initial gravitational field. As time goes on, all of these displacements grow by the same overall factor, the linear growth factor D(t). The whole evolution is captured in one equation:

x(q, t) = q + D(t) · s(q)

Here x is the particle's comoving position now, q is where it started, D(t) is the time-dependent growth factor (D = 0 at early times, growing monotonically), and s(q) is the time-independent displacement field, fixed once and for all by the initial conditions through s(q) = −∇_qΦ for an initial potential Φ. Every particle knows from birth where it is heading; gravity only sets the speed of the parade, never the direction.

That single property — a frozen displacement direction — is what makes pancakes inevitable, and it is the heart of the whole scheme.

How it works: the deformation tensor and the order of collapse

To see why sheets come first, linearize the displacement around a point. The local stretching and squeezing of the matter is governed by the deformation tensor, the gradient of the displacement field:

D_ij = ∂s_i / ∂q_j

This is a symmetric 3×3 tensor (for an irrotational, potential flow), so at each point it has three real eigenvalues, conventionally ordered

λ1 ≥ λ2 ≥ λ3

with eigenvectors that define the local principal axes. In the eigenframe the Zel'dovich map factorizes, and the local density follows directly from conservation of mass — the new density is the initial density divided by the Jacobian of the map:

ρ(x, t) = ρ̄ / [ (1 − D(t) λ1)(1 − D(t) λ2)(1 − D(t) λ3) ]

Now read off the physics. Each factor (1 − D λ_i) shrinks toward zero as D grows, and the density blows up when any factor reaches zero. The first factor to vanish is the one with the largest positive eigenvalue, λ1, at the moment D(t) = 1/λ1. At that instant the matter has collapsed to zero thickness along a single axis — the λ1 eigendirection — while the other two directions are still extended. The result is a two-dimensional sheet: Zel'dovich's pancake.

Only later, as D climbs further to 1/λ2, does the sheet collapse along its second axis to make a filament; later still, at D = 1/λ3, the filament pinches into a knot. The order is forced: pancakes form first, then filaments, then nodes. A perfectly spherical collapse (λ1 = λ2 = λ3) would skip the sheet stage, but that requires three eigenvalues to be exactly equal — a measure-zero coincidence that essentially never happens in a random field. Anisotropy is generic, so sheets are generic.

A worked example: when does a patch pancake?

Consider a region of the early universe whose deformation tensor, smoothed on some scale, has eigenvalues

λ1 = 1.5,   λ2 = 0.6,   λ3 = −0.4   (in units where the field is normalized to today, D(t_0) = 1)

Pancake collapse along the λ1 axis happens when D(t) = 1/λ1 = 1/1.5 ≈ 0.67. In an Einstein–de Sitter universe the growth factor scales with the scale factor, D ∝ a, so D = 0.67 corresponds to a ≈ 0.67, i.e. redshift z = 1/a − 1 ≈ 0.5. This patch sheeted up when the universe was about two-thirds its present size.

The second axis collapses at D = 1/λ2 = 1/0.6 ≈ 1.67 — which is larger than 1, the value today. So this patch is a sheet now and will only thread into a filament in the future. The third eigenvalue is negative (λ3 = −0.4), so that factor, (1 − D λ3) = (1 + 0.4 D), keeps growing; that axis is expanding, never collapsing. This region is a wall feeding a future filament, with one direction permanently draining outward into a void. That eigenvalue bookkeeping — count how many of the three (1 − D λ_i) factors have gone negative — is exactly how modern cosmic-web classifiers (the "T-web") label a point as void (0 collapsed axes), sheet (1), filament (2), or knot (3).

The headline numbers are worth keeping: the recipe is x = q + D(t)·s(q); the order of collapse always puts pancakes first; and the entire construction is exact to first order and valid right up to shell crossing — the moment the first factor (1 − D λ1) tries to go negative.

Shell crossing: the wall the approximation hits

The density formula above does something alarming at D = 1/λ1: it diverges, and immediately afterward the factor (1 − D λ1) goes negative, making the density negative. That is not physics; it is the signal that the approximation has reached its domain boundary. The event is shell crossing: the map q → x has stopped being one-to-one. Two mass elements that launched from different starting positions q now arrive at the same physical point x, and the Jacobian determinant of the map passes through zero.

What really happens at shell crossing is that particle trajectories intersect. In a true gravitating system, the streams would interpenetrate and the region would become "multi-stream" — at a given point you would find matter moving in several directions at once, oscillating back and forth across the dense sheet, building a stable caustic. The Zel'dovich approximation does not know about gravity after this point; it keeps every particle moving in a straight line at constant comoving direction. So the freshly assembled pancake immediately gets smeared out as particles fly straight through and out the far side. The sheet you worked so hard to build evaporates.

This is the single most important caveat: the Zel'dovich approximation is a beautiful, accurate description of how structure assembles right up to the first crossing, and an unphysical one after. Everything that follows in the field — 2LPT, the adhesion model, full N-body — is a way of either delaying that wall or replacing the post-crossing behavior with something physical.

Variants and extensions

Scheme	What it adds	Order	Behavior at/after shell crossing	Typical use
Eulerian linear theory	δ = D(t)·δ₀ at fixed points	1st (Eulerian)	Never collapses; density stays small	Large-scale δ, power spectrum growth
Zel'dovich (ZA / 1LPT)	Displacement x = q + D·s(q)	1st (Lagrangian)	Exact until crossing, then smears pancakes	N-body ICs, intuition, web classification
2LPT	Second-order displacement term ∝ D²	2nd (Lagrangian)	Delays crossing; better mildly nonlinear field	Standard modern initial conditions
3LPT / ALPT	Third order / augmented hybrid	3rd	Further delay; spheroidal collapse hybrids	Fast mock catalogs (PATCHY, etc.)
Adhesion model	Artificial viscosity (Burgers' equation)	ZA + viscosity	Traps matter at caustics; stable web skeleton	Cosmic-web morphology, void/filament maps
Full N-body	Solves gravity for all particles	Nonlinear (exact)	Correct multi-stream caustics, virialization	Precision predictions, halo catalogs

The progression is a ladder of fidelity. Eulerian linear theory is below Zel'dovich because it cannot collapse anything at all. Zel'dovich captures the first collapse exactly. 2LPT — second-order Lagrangian perturbation theory — adds the next term in the displacement expansion (a piece growing like D²) and pushes shell crossing to later times, which is why almost every modern N-body simulation starts from 2LPT rather than bare Zel'dovich. The adhesion model takes a completely different tack: it modifies the equation of motion with an artificial viscosity so that the Zel'dovich flow becomes the inviscid limit of Burgers' equation, whose shocks trap matter at caustics instead of letting it stream through. The result is a sharp, stable cosmic-web skeleton that matches simulations remarkably well without any free-streaming smear.

Common pitfalls and misconceptions

• Confusing the Lagrangian and Eulerian pictures. Eulerian linear theory writes δ(x) = D(t)·δ₀(x): the density at a fixed point grows uniformly and nothing ever collapses. Zel'dovich is Lagrangian — it moves the particles. Tracking trajectories, not fields, is precisely what lets density become singular and pancakes appear.
• Thinking collapse is spherical. The "spherical collapse model" is a useful toy, but generic patches are triaxial. Because one eigenvalue is almost always strictly largest, collapse is anisotropic and sheet-first. Spherical collapse is the exception, not the rule.
• Trusting the density formula past shell crossing. Once (1 − D λ1) hits zero the formula gives infinite and then negative densities. That is a hard boundary, not a feature to extrapolate through. The approximation is silent about virialized halos.
• Assuming the displacement direction evolves. It does not — s(q) is fixed by the initial conditions; only the scalar D(t) changes. The straight-line motion is the source of both the scheme's elegance and its post-crossing failure.
• Equating "pancake" with a real, persistent object. In bare ZA the pancake is transient: it forms then smears. Persistent walls require the adhesion model or full gravity, which trap matter at the caustic.
• Forgetting it is a first-order theory. The accuracy (≈10–20% in the power spectrum on mildly nonlinear scales) is good but not exact; for precision work you climb the LPT ladder or run N-body.

Observational status and applications

Zel'dovich's pancake prediction was a genuine forecast: in 1970 the large-scale distribution of galaxies was essentially unmapped, and a sheet-and-void universe was far from the default expectation. The vindication came from redshift surveys. The 1986 CfA survey's "Stick Man" slice revealed walls and voids; the Las Campanas, 2dF, and Sloan Digital Sky Surveys mapped the full cosmic web of sheets, filaments, and clusters threaded by enormous voids — exactly the cellular morphology Zel'dovich described.

The approximation also remains a working tool, not just a historical milestone:

• Initial conditions for simulations. Every cosmological N-body run begins by displacing a uniform particle grid with a Zel'dovich (or 2LPT) displacement field drawn from the input power spectrum. Codes like 2LPTic, MUSIC, and monofonIC do exactly this.
• BAO reconstruction. Galaxy surveys (BOSS, eBOSS, DESI) sharpen the baryon acoustic oscillation standard ruler by "reconstructing" — moving galaxies backward along the negative Zel'dovich displacement to partially undo nonlinear smearing. This roughly halves the error on the acoustic scale and hence on the distance–redshift relation, directly improving dark-energy constraints.
• Fast mock catalogs. Approximate methods (PATCHY, COLA, FastPM) lean on Lagrangian perturbation theory to generate thousands of mock universes cheaply for covariance estimation, where full N-body would be prohibitively expensive.
• Cosmic-web classification. The eigenvalue-counting logic of the deformation tensor underlies the T-web, V-web, and NEXUS schemes that segment the universe into voids, sheets, filaments, and knots.

Quantitative analysis: from continuity to the pancake

Start from mass conservation. The mass in a Lagrangian volume element d³q must equal the mass in the physical element d³x it maps to:

ρ(x) d³x = ρ̄ d³q   ⟹   ρ(x) = ρ̄ / |det(∂x_i/∂q_j)|

The Jacobian of the Zel'dovich map x = q + D·s(q) is

∂x_i/∂q_j = δ_ij + D(t) (∂s_i/∂q_j) = δ_ij + D(t) D_ij

Diagonalizing the deformation tensor D_ij into its eigenvalues λ1, λ2, λ3, the determinant factorizes and we recover the density formula:

ρ(x, t) = ρ̄ / [ (1 + D λ̃1)(1 + D λ̃2)(1 + D λ̃3) ]

where the λ̃ are the eigenvalues of ∂s/∂q (with the sign convention chosen so a positive eigenvalue means convergence — collapse — which flips the sign to the (1 − D λ) form used above). In the linear, small-displacement limit the density contrast is δ = ρ/ρ̄ − 1 ≈ −D(λ̃1 + λ̃2 + λ̃3) = −D ∇·s, recovering Eulerian linear growth δ ∝ D as it must — but the full factored form keeps going where the linear expansion cannot, all the way to the singularity.

The velocity follows by differentiating the position with respect to time. Defining the conformal-time growth, the peculiar velocity is

v(q, t) = a Ḋ(t) s(q)   ∝   ḋ · displacement

so the velocity field is everywhere parallel to the (fixed) displacement field, with a magnitude set by the growth rate. This is the defining feature: velocity and displacement point the same way and never rotate. The Zel'dovich flow is curl-free (irrotational) by construction, which is why a single scalar potential s = −∇Φ suffices to describe the entire kinematics — until trajectories cross and the smooth potential description breaks down.