Solving the Navier-Stokes Equations: What "Solution" Really Means

When someone asks whether the Navier-Stokes equations can be solved, the honest answer is: it depends on what you mean by "solved." Engineers solve them numerically every day to design aircraft and predict weather. Mathematicians have known exact pencil-and-paper solutions for simple geometries since the 1840s. But nobody has solved the problem, which is proving that solutions always exist and behave well.

There are three distinct levels:

• Exact solutions. In special geometries, symmetry simplifies the equations so much that you can write the answer down in closed form.
• Numerical solutions. For real-world geometries, computers approximate the flow cell by cell and time step by time step. This is computational fluid dynamics (CFD).
• The existence question. Does every smooth starting flow in 3D stay smooth forever? That's the Millennium Prize problem, and it's open.

The rest of this page walks through each level and, more importantly, the gaps between them. The gaps are where the interesting mathematics lives.

Fix the incompressible system $$\partial_t u + (u\cdot\nabla)u = -\frac{1}{\rho}\nabla p + \nu\,\Delta u + f, \qquad \nabla\cdot u = 0.$$ "Solving" it can mean three inequivalent things:

1. Closed form: exhibiting fields $u(x,t)$ and $p(x,t)$ that satisfy the system exactly. In practice this requires a symmetry-driven reduction that turns the PDE system into a linear ODE or a heat-type equation.
2. Discrete approximation: constructing $u_h$ on a mesh of size $h$ such that $u_h$ converges to a solution as $h\to 0$, subject to stability and consistency constraints on the scheme.
3. Well-posedness: the Clay Mathematics Institute formulation asks whether, for every smooth, divergence-free, rapidly decaying initial datum $u_0$ on $\mathbb{R}^3$ (or smooth data on the torus), there is a smooth solution $(u,p)$ with bounded energy for all $t>0$.

Each level presupposes less structure than the previous one and asks for more generality. The first two are settled practice. The third is the open problem.

Sometimes the geometry is so simple that the hardest part of the equations disappears. In a straight pipe or between two parallel plates, the fluid moves in one direction and varies in another. That layout kills the nonlinear term entirely, and what's left is an equation you can integrate by hand.

The classical families every fluids course covers:

• Poiseuille flow: pressure-driven flow through a pipe, with its parabolic velocity profile and fourth-power flow law
• Couette flow: fluid sheared between a moving plate and a fixed one
• Stokes's first problem: fluid set in motion by a suddenly moving wall
• The Lamb-Oseen vortex: a single vortex slowly spreading out as viscosity acts
• Stagnation-point flow: flow hitting a wall head-on

These aren't museum pieces. They're benchmarks for numerical methods, and they carry real intuition about how viscosity shapes a flow. Our exact solutions page compares the main families, while the Poiseuille guide gives the full pipe-flow derivation.

For a parallel flow $u=(u(y,t),0,0)$ the divergence-free condition holds automatically and the advection term vanishes identically: $(u\cdot\nabla)u = 0$. The momentum equation collapses to $$\frac{\partial u}{\partial t} = -\frac{1}{\rho}\frac{dp}{dx} + \nu\,\frac{\partial^2 u}{\partial y^2},$$ a linear equation. In steady state it becomes a second-order ODE solvable by direct integration; between plates at $y=0$ and $y=h$ this yields the plane Poiseuille profile $$u(y)=\frac{1}{2\mu}\left(-\frac{dp}{dx}\right)y(h-y).$$

The circular-pipe reduction uses the radial operator $r^{-1}(ru')'$ instead. The complete Hagen–Poiseuille derivation gives $u(r)=G(R^2-r^2)/(4\mu)$ and $Q=\pi GR^4/(8\mu)$. The pattern generalizes: every classical exact solution exists because some symmetry (parallel streamlines, radial symmetry, self-similarity) removes or tames the nonlinear feedback between the velocity field and its own transport. Break the symmetry and the closed form is gone. That is why exact solutions, valuable as they are, say almost nothing about the general three-dimensional problem.

For anything more complicated than a straight pipe (a wing, a blood vessel, a storm system) there's no formula. Instead, CFD chops space into millions of small cells and marches the flow forward in time, computing how each cell's velocity and pressure respond to its neighbors. You don't solve the original continuous equations. You solve a discrete stand-in for them, built so that its answers approach the true flow as the cells shrink.

The main discretization families:

• Finite difference: replace derivatives with differences on a grid
• Finite volume: track fluxes through cell faces so conserved quantities stay conserved (most commercial CFD works this way)
• Finite element: build the solution from simple basis functions
• Spectral methods: expand the flow in Fourier modes, extremely accurate for periodic domains
• Lattice Boltzmann: simulate simplified particle populations whose averages recover the fluid equations

Turbulent flows add another layer of modeling choices, from resolving everything (DNS) to averaging everything (RANS). How fine you have to resolve depends on the Reynolds number, and which equations you discretize depends on whether the flow is incompressible or compressible.

The central structural difficulty is pressure-velocity coupling. Incompressible Navier-Stokes has no evolution equation for $p$; pressure acts as a Lagrange multiplier enforcing $\nabla\cdot u = 0$. Standard treatments include projection methods (Chorin, 1968), which advance a provisional velocity and then project onto divergence-free fields via a pressure Poisson solve, pressure-correction schemes of the SIMPLE/SIMPLEC/PISO family, and fully coupled solvers.

Discretizing advection forces a tradeoff: upwind schemes are stable but introduce numerical diffusion, while central schemes are more accurate but prone to spurious oscillations. Time integration divides into explicit schemes (cheap per step, CFL-limited), implicit schemes (stable but requiring nonlinear solves per step), and IMEX combinations.

For turbulence the hierarchy is DNS, LES, RANS. Direct numerical simulation resolves all dynamically active scales, with cost growing so fast in Reynolds number that it's confined to modest $Re$. Large-eddy simulation resolves large scales and models subgrid stresses. Reynolds-averaged formulations model the entire fluctuation field and close the resulting moment equations with empirical models.

Here's the twist. Engineers compute fluid flows daily and physicists trust the equations completely, yet mathematicians can't prove the equations always give a well-behaved answer. Start from any smooth flow pattern in 3D and let it evolve. Does it stay smooth forever, or can it concentrate energy so violently that velocities become infinite in finite time? Nobody knows. That question, existence and smoothness, is one of the seven Millennium Prize Problems, worth $1 million.

What's already settled:

• In 2D, the answer is yes: smooth solutions exist for all time. That was worked out in the mid-20th century.
• In 3D for short times, smooth solutions exist. The unknown is whether "short" can always be extended to "forever."
• In 3D for all time, a weaker kind of solution (Leray's weak solutions) always exists. Whether those weak solutions are actually smooth is exactly the open question.

Why 3D resists when 2D doesn't comes down to vortex stretching and the way energy moves to small scales; our page on why the problem is hard goes deeper.

Leray (1934) proved global-in-time existence of weak solutions for divergence-free initial data $u_0\in L^2(\mathbb{R}^3)$: fields satisfying the equations in the distributional sense and the energy inequality, but with no guaranteed pointwise regularity. Uniqueness of these solutions is unknown.

The Clay problem asks: for $f\equiv 0$ and smooth, divergence-free, rapidly decaying $u_0$ on $\mathbb{R}^3$, does there exist a pair $(u,p)$ with $u\in C^\infty(\mathbb{R}^3\times[0,\infty))$ and globally bounded energy? Either a proof of global regularity or a demonstration of finite-time breakdown for admissible data resolves the problem.

The dimensional split is sharp. Global regularity in 2D was established through the combined work of Leray, Ladyzhenskaya, Lions, and Prodi; the key is that 2D vorticity is transported without stretching, giving strong a priori control. In 3D, local-in-time existence of smooth solutions is classical, and smoothness is known to persist as long as certain scale-critical norms stay finite, but global regularity for general data remains open in both directions: no proof, no counterexample.

If computers can simulate any flow you like, why does the mathematical question matter? Because a simulation and a proof answer different questions. A simulation gives you numbers on a grid. A proof gives you a guarantee about the equations themselves: that solutions exist, are unique, and respond continuously to changes in the starting data. Without that guarantee there's a logical gap at the foundation of fluid mechanics. "The equations work" and "we can prove the equations work" aren't the same statement.

The gap has practical teeth. Simulations sometimes blow up, and when one does, you can't always tell whether you've hit a numerical artifact (mesh too coarse, time step too large) or brushed against a genuine singularity of the equations themselves. A regularity theory would settle which failures are real. The strategies aimed at closing this gap live on our approaches and subproblems pages.

Convergence theorems for numerical schemes are conditional: they typically assume the limiting solution has exactly the regularity the Millennium problem asks us to establish. Absent that, mesh refinement produces sequences whose limits are guaranteed to be weak solutions at best, and the weak class has surprises. Convex-integration methods (developed for Euler in the resolution of Onsager's conjecture by Isett and by Buckmaster, De Lellis, Székelyhidi, and Vicol, and extended by Buckmaster and Vicol to Navier-Stokes in classes weaker than Leray-Hopf) show that distributional solutions can be badly non-unique. The notion of "the" computed solution is therefore less innocent than it looks.

There's also a resolution gap. DNS currently reaches $Re\sim 10^4\text{-}10^5$ on large machines, while atmospheric flows operate at $Re\sim 10^9$. Engineering practice extrapolates trust in the equations across several orders of magnitude in scale separation, without a rigorous well-posedness theory underwriting the extrapolation. None of this makes CFD unreliable in practice; it locates precisely what a proof would add.

Where things stand today:

• Exact solutions: a settled catalog for symmetric flows, still growing slowly at the margins
• Numerical solutions: mature and spectacularly accurate for engineering, with reliability that degrades as turbulence modeling carries more of the load
• The existence question: open. No proof of global smoothness, no blowup counterexample, prize unclaimed

The hierarchy of "solving" mirrors the hierarchy of understanding. We can compute nearly any flow and approximate it as finely as budgets allow, but we can't yet prove the most basic well-posedness statement about the 3D equations. Partial results keep arriving, and our progress page tracks them as they land.

The strongest unconditional 3D result remains Caffarelli-Kohn-Nirenberg partial regularity: for suitable weak solutions, the singular set has one-dimensional parabolic Hausdorff measure zero. Conditional regularity criteria (Beale-Kato-Majda via vorticity, originally formulated for Euler and adapted to Navier-Stokes, and the Prodi-Serrin-Ladyzhenskaya scale-critical integrability conditions) convert control of a single critical quantity into full smoothness.

Active directions include convex integration and the structure of non-uniqueness below the Leray-Hopf class, quantitative regularity and blowup-rate results, computer-assisted analysis of candidate self-similar blowup scenarios, and stochastic regularization questions. The scaling gap (the known a priori estimates are supercritical relative to the natural scaling of the equations) is still the central obstruction, and it's why incremental improvements haven't closed the problem.