Euler Fluid Equations

The Euler fluid equations are a coupled system of quasilinear first-order partial differential equations that govern the motion of an idealized inviscid and adiabatic continuum fluid—one in which molecular viscosity (μ = 0) and thermal conductivity (κ_th = 0) vanish. Published by Swiss mathematician and physicist Leonhard Euler in 1757 in his Berlin Academy memoir Principes généraux du mouvement des fluides, the equations represent the earliest complete partial differential formulation of continuum mechanics, predating the viscous Navier-Stokes equations by nearly seven decades.

Mathematically, the Euler equations emerge as the formal high-Reynolds number limit (Re → ∞) of the Navier-Stokes system away from solid walls and thin dissipation sheets governed by boundary layer theory. Because dropping second-order viscous diffusion transforms the governing equations into a first-order hyperbolic system of conservation laws, the Euler equations capture the propagation of acoustic waves, gravity and Rossby waves in geophysical fluid dynamics, steady aerodynamic potential flow governed by Bernoulli's principle, and non-linear compressible wave steepening into supersonic shock waves across aerospace engineering and astrophysics.

Differential and Conservative Forms of the Euler System

In an Eulerian coordinate frame x = (x, y, z), the state of a compressible single-component fluid is specified by five primitive field variables: mass density ρ(x, t), three velocity components u(x, t) = (u, v, w), and thermodynamic pressure p(x, t) (closed via an equation of state with internal energy per unit mass e(ρ, p)). Applying conservation of mass, Newton's second law (with isotropic normal Cauchy stress σ = −pI and external body acceleration g), and the first law of thermodynamics to an infinitesimal fluid parcel yields Euler's three governing equations:

• Conservation of Mass (Continuity Equation):

• Conservation of Inviscid Momentum:

• Conservation of Total Energy:

where E = ρ(e + ½|u|²) is the total energy per unit volume (internal plus kinetic) and H₀ = (E + p)/ρ = h + ½|u|² is the stagnation enthalpy per unit mass. For an ideal calorically perfect gas with specific heat ratio (adiabatic index) γ = c_p/c_v (where γ = 1.4 for diatomic air), the equation of state closes the system via p = (γ − 1)ρe = ρRT. Where the solution remains smooth and shock-free, the energy equation is equivalent to Lagrangian conservation of specific entropy s along particle paths (isentropic flow: Ds/Dt = ∂s/∂t + u · ∇s = 0, or D(pρ^−γ)/Dt = 0).

In divergence (conservation) form—required when integrating across discontinuous shock waves—the system is written compactly as ∂U/∂t + ∇ · F(U) = S, where the conserved state vector U and flux tensor F are:

Singular Perturbation Limit of the Navier-Stokes Equations

Although the Euler momentum equation is obtained formally by setting dynamic viscosity μ = 0 (or 1/Re = 0) in the Navier-Stokes equations, this reduction is a singular perturbation because it eliminates the second-order Laplacian operator (μ∇²u), lowering the spatial differential order of the partial differential equation from second order to first order:

Irrotational Flow and Bernoulli's Theorems

Using the vector identity (u · ∇)u = ∇(½|u|²) + ω × u (where ω = ∇ × u is vorticity), the Euler momentum equation in a conservative gravitational field (g = −∇Φ) takes Crocco's form:

By Kelvin's circulation theorem and the vorticity equation, an inviscid barotropic flow that starts from a uniform or irrotational state (ω = ∇ × u = 0) remains irrotational for all time provided no shock waves form. Because ∇ × u = 0, the velocity field can be expressed as the gradient of a scalar velocity potential φ(x, t): u = ∇φ. Substituting u = ∇φ into the Crocco–Euler equation and integrating spatially yields the unsteady Bernoulli equation valid throughout the entire fluid domain:

In steady flow (∂/∂t = 0), this reduces to Bernoulli's principle along streamlines (or everywhere in irrotational flow): ½u² + p/ρ + gz = constant (for incompressible flow) or ½u² + γ/(γ − 1) + gz = constant (for compressible isentropic flow).

D'Alembert's 1752 Paradox of Zero Drag

When an incompressible (∇ · u = 0) AND irrotational (u = ∇φ) fluid flows past a submerged body, substituting u = ∇φ into the continuity equation yields Laplace's equation:

subject to the inviscid no-penetration boundary condition ∂φ/∂n = u · n = 0 on the body surface and ∇φ → U_∞ at infinity. In 1752, French mathematician Jean le Rond d'Alembert proved that for any steady, unseparated potential flow around a closed finite body in an unbounded fluid, the fore-and-aft surface pressure distribution computed from Bernoulli's equation is symmetric: high stagnation pressure on the nose is balanced by an equal pressure recovery at the rear stagnation point, resulting in identically zero net hydrodynamic force (F = 0)—zero drag and zero lift.

This stark contradiction with everyday experience—known as d'Alembert's paradox—persisted for 152 years until Ludwig Prandtl's 1904 boundary layer theory demonstrated that no matter how small the viscosity ν, the no-slip boundary condition creates a thin wall boundary layer that separates under the rear adverse pressure gradient, forming a low-pressure rotational wake that destroys fore–aft pressure symmetry.

Linear Acoustic Waves and the Speed of Sound

To analyze small-amplitude pressure disturbances in a compressible fluid at rest (ρ₀, p₀, u₀ = 0), one expands the state variables as ρ = ρ₀ + ρ', p = p₀ + p', and u = u' with |ρ'| ≪ ρ₀. Neglecting quadratic products of primed perturbations in the Euler continuity and momentum equations yields the linearized acoustic system:

Because acoustic compressions occur on time scales far faster than thermal conduction, the process is isentropic, so pressure and density fluctuations are linked by p' = (∂p/∂ρ)_s ρ' = c_s² ρ'. Taking the time derivative of the continuity equation and substituting the divergence of the momentum equation produces the linear wave equation:

where the phase speed of acoustic waves (first derived correctly by Pierre-Simon Laplace in 1816 by accounting for adiabatic heating over Isaac Newton's isothermal formula) is the adiabatic speed of sound:

Wave Steepening, Rankine–Hugoniot Jump Conditions, and Riemann Solvers

For finite-amplitude compressible flows (where the Mach number Ma = |u|/c_s is no longer infinitesimal), the Euler equations in one spatial dimension admit three real characteristic propagation speeds (eigenvalues of the flux Jacobian ∂F/∂U):

Along a forward-propagating Riemann wave, both local particle velocity u and local sound speed c_s (which scales with T^1/2 and ρ^(γ−1)/2) are higher at the compressed wave crest than in the trough ahead of it. Consequently, wave crests travel faster (u + c_s) than wave troughs, causing the waveform to steepen progressively until ∂u/∂x → −∞ in finite time (the "gradient catastrophe" analyzed by Bernhard Riemann in 1860), forming a discontinuous shock wave.

Because spatial derivatives do not exist at a discontinuity, weak solutions of the Euler equations across a shock moving at speed S_s are governed by the integral conservation form ∂U/∂t + ∂F/∂x = 0, which yields the algebraic Rankine–Hugoniot jump conditions derived by William John Macquorn Rankine (1870) and Pierre-Henri Hugoniot (1887). In a reference frame stationary with respect to a normal shock (linking upstream supersonic state 1 with Ma₁ > 1 to downstream subsonic state 2 with Ma₂ < 1):

• Mass Flux Continuity: ρ₁ u₁ = ρ₂ u₂
• Momentum Flux Balance: p₁ + ρ₁ u₁² = p₂ + ρ₂ u₂²
• Stagnation Enthalpy Continuity: h₁ + ½u₁² = h₂ + ½u₂²
• Density & Pressure Jump Ratios: ρ₂/ρ₁ = u₁/u₂ = [(γ + 1)Ma₁²] / [2 + (γ − 1)Ma₁²] (bounded above by (γ + 1)/(γ − 1) = 6 for air as Ma₁ → ∞), and p₂/p₁ = 1 + 2γ/(γ + 1).

To prevent unphysical expansion shocks that violate the second law of thermodynamics, weak Euler solutions must satisfy the Lax entropy condition (s₂ > s₁, requiring Ma₁ ≥ 1). In modern computational fluid dynamics (CFD) and astrophysical hydrodynamics, finite-volume Godunov methods—such as the Roe, HLLC, and exact Riemann solvers—solve the local piecewise-constant Euler jump problem across every grid cell interface at each time step to capture shocks without spurious numerical oscillations.

How can the inviscid Euler equations produce irreversible entropy jumps across shock waves if viscosity and heat conduction are zero?

For smooth differentiable velocity fields, the Euler equations strictly conserve entropy (Ds/Dt = 0). Once nonlinear wave steepening creates a shock discontinuity, however, spatial derivatives diverge and the differential equivalence between total energy conservation and entropy conservation breaks down. Enforcing integral conservation of mass, momentum, and total energy across the jump (the Rankine–Hugoniot conditions) forces kinetic energy to convert into internal thermal energy across the shock front, producing a finite entropy increase Δs ∝ (Ma₁² − 1)³ for weak shocks. Physically, this jump represents the vanishing-viscosity limit (μ, κ_th → 0) of a real Navier-Stokes shock layer whose physical thickness δ_shock ∼ ν/c_s shrinks to zero while its integrated dissipation remains independent of μ (analogous to the dissipative anomaly in the Kolmogorov turbulence cascade).

Why must the Euler equations be written in conservative (divergence) form when simulating supersonic flows?

In regions where flow variables are smooth, the non-conservative primitive form (∂u/∂t + u · ∇u = −∇p/ρ) and the conservative divergence form (∂(ρu)/∂t + ∇ · (ρu ⊗ u + pI) = 0) are mathematically identical via the product rule. Across a shock discontinuity, the product rule fails; by the Lax–Wendroff theorem (1960), numerical discretizations based on the non-conservative primitive equations converge to weak solutions with incorrect shock propagation speeds and violated Rankine–Hugoniot jump balances.

How are the Euler equations used in weather forecasting and geophysical fluid dynamics?

Outside the lowest ~1 km planetary boundary layer, molecular viscosity has negligible direct impact on synoptic-scale atmospheric motions (Re ∼ 10¹²). Consequently, the non-hydrostatic dynamical cores of modern numerical weather prediction models (such as the Weather Research and Forecasting model, WRF, and MPAS) solve the three-dimensional compressible Euler equations on a rotating sphere—augmented by the Coriolis effect, moisture phase-change thermodynamics, and parameterized sub-grid turbulent eddy diffusion.

1. Euler, L. (1757). "Principes généraux du mouvement des fluides." Mémoires de l'Académie des Sciences de Berlin, 11, 274–315. English translation by T. E. Burton & U. Frisch in Physica D: Nonlinear Phenomena, 237(14–17), 1825–1839 (2008). doi:10.1016/j.physd.2008.04.015
2. Lax, P. D. (1973). Hyperbolic Systems of Conservation Laws and the Mathematical Theory of Shock Waves. Regional Conference Series in Applied Mathematics, Vol. 11, SIAM. doi:10.1137/1.9781611970562
3. Toro, E. F. (2009). Riemann Solvers and Numerical Methods for Fluid Dynamics: A Practical Introduction (3rd ed.). Springer-Verlag. doi:10.1007/b79761
4. Majda, A. J., & Bertozzi, A. L. (2002). Vorticity and Incompressible Flow. Cambridge Texts in Applied Mathematics, Cambridge University Press. doi:10.1017/CBO9780511613203