Smoothed-Particle Hydrodynamics

In 1977, two independent teams — Gingold & Monaghan and Lucy — needed to simulate astrophysical gas clouds without a fixed computational grid. Grids struggle when matter tears apart or collapses; the teams' solution was radical: make the fluid itself carry the equations. The result, Smoothed-Particle Hydrodynamics (SPH), is still in active use half a century later.

The core idea is deceptively simple. Instead of carving space into cells, you scatter particles through the fluid. Each particle carries its own mass, velocity, density and pressure. To evaluate any quantity at a point in space, you sum the contributions of nearby particles, each weighted by a smoothing kernel — a bell-shaped function $W(r, h)$ that falls to zero beyond a radius $h$ called the smoothing length.

Because the particles move with the fluid, SPH is Lagrangian: it tracks fluid elements rather than fixed points in space. There is no mesh to tangle, no cells to refill — just particles following their own trajectories while continuously averaging over their neighbors.

This article explains how that averaging works, what it costs, and why SPH shows up in everything from numerical simulations of colliding stars to the water in your favorite blockbuster film.

The demo below runs a miniature SPH simulation entirely in your browser. A blob of particles drops under gravity, hits the floor, and settles. Each particle pushes its neighbors through a pressure force proportional to the density gradient computed by the smoothing kernel.

<!-- {{c_html_comment}} -->
<div class="controls">
  <label>{{lbl_gravity}} <input type="range" id="gravity" min="100" max="800" value="400" step="10"></label>
  <label>{{lbl_damping}} <input type="range" id="damping" min="0" max="99" value="20" step="1"></label>
  <button id="btn-drop" type="button">{{btn_drop}}</button>
  <button id="btn-reset" type="button" class="ghost">{{btn_reset}}</button>
</div>
<canvas id="c" width="480" height="340"></canvas>
<div id="status" class="status">{{status_idle}}</div>

/* {{c_css_comment}} */
* { box-sizing: border-box; }
body { font-family: system-ui, sans-serif; margin: 0; background: #f0f4f8; }
.controls { display: flex; gap: .6rem; flex-wrap: wrap; align-items: center;
            padding: .5rem .4rem; background: #e2eaf2; border-radius: 8px; margin-bottom: .5rem; }
label { font-size: .82rem; color: #334; display: flex; align-items: center; gap: .3rem; }
input[type=range] { width: 80px; cursor: pointer; }
canvas { display: block; border-radius: 8px; border: 1px solid #c4d0dc; background: #dce8f5; width: 100%; }
.status { font-size: .88rem; font-weight: 600; color: #2a4d6e; margin-top: .4rem; min-height: 1.3em; }
button { font: 600 13px system-ui; padding: .38rem .8rem; border: 1px solid #2a4d6e;
         background: #2a4d6e; color: #fff; border-radius: 6px; cursor: pointer; }
button.ghost { background: #fff; color: #2a4d6e; }

// Code not found

Try adjusting gravity and damping with the sliders. Notice how increasing the smoothing radius makes the fluid look more viscous — particles "feel" more neighbors and average out sudden pressure spikes. Decrease it too much and the simulation becomes noisy: each particle sees too few neighbors and pressure estimates jump wildly. The sweet spot is roughly 30–50 neighbors per particle, a rule of thumb SPH practitioners have used since the 1980s.

SPH is elegant, but naive implementations hide a cost that scales badly.

One time step breaks into three phases:

1. Neighbor search — for each particle, find every other particle within radius $h$.
2. Density estimation — sum $\rho_i = \sum_j m_j W(|\mathbf{r}_i - \mathbf{r}_j|, h)$ over neighbors.
3. Force evaluation — apply the SPH form of the Navier-Stokes equations: pressure gradient, viscosity and any external forces.

With $N$ particles, the naive pairwise search costs $O(N^2)$ per step — fine for a few hundred particles, disastrous at a million.

The fix: spatial hashing. Divide space into cells of size $h$. Each particle hashes to its cell; to find neighbors you only look at the $3^d$ adjacent cells ($3^2 = 9$ in 2-D, $3^3 = 27$ in 3-D). If particles are roughly uniformly spread, the neighbor count per particle is bounded by a constant, dropping the per-step cost to $O(N)$ for force evaluation. The hash table rebuild costs $O(N)$ with a counting sort, giving an overall per-step complexity of $O(N \log N)$ in practice (dominated by sort on GPU) or $O(N)$ in ideal uniform distributions.

Accuracy scales with both $N$ (more particles) and $h$ (smaller kernel). Reducing $h$ improves spatial resolution but forces smaller time steps (the CFL condition: $\Delta t \leq C \cdot h / v_{\max}$), so the total cost to reach time $T$ scales roughly as $O(N^{4/3})$ in three dimensions for fixed accuracy — an uncomfortable truth for high-resolution ocean simulations.

SPH sits in a crowded zoo of meshfree methods alongside the Material Point Method (MPM) and Moving Least Squares approaches. Each trades accuracy, stability and cost differently, but the neighbor-search bottleneck is universal: no matter how you smooth, you must still find who is close to whom.

Because SPH needs no grid, it excels wherever geometry changes violently:

• Astrophysics: the method's birthplace. SPH simulates galaxy mergers, proto-planetary disk formation and supernova shockwaves. The GADGET and AREPO codes run SPH at cosmological scales.
• Tsunami and dam-break modeling: water that tears across irregular terrain is exactly the free-surface, large-deformation regime where grids fail and SPH thrives.
• Metal forming and impact mechanics: solid materials undergoing fracture and plastic flow are treated as nearly-incompressible fluids in SPH, modeling bullet impacts, explosions and crash tests.
• Visual effects: DreamWorks, ILM and Sony Pictures Imageworks have used SPH-based solvers (often coupled with FLIP) to animate water, lava and crowds of tumbling objects in films and games.
• GPU acceleration: SPH maps cleanly onto massively parallel hardware. NVIDIA's PhysX and Bullet Physics both include GPU SPH for real-time fluid in games.

The same kernel-interpolation idea that makes SPH work for fluids also appears in kernel density estimation in statistics and in radial basis function interpolation — a reminder that the mathematical core of SPH is more general than any single application domain. It connects naturally to the broader landscape of numerical simulations and optimization problems that drive scientific computing.

SPH's bet is simple: if the fluid is the solver, there is no mesh to corrupt, no cell to flood, no topology to repair. Particles carry mass, momentum and the smoothing kernel does the rest — integrating the Navier-Stokes equations one neighbor-weighted sum at a time.

That bet has paid off across fifty years and a dozen disciplines. From the galaxies Gingold and Monaghan wanted to collide in 1977, to the oceans of the latest blockbuster, SPH's meshfree Lagrangian approach keeps proving itself wherever matter refuses to stay still.

The price is the neighbor search: as $N$ grows, keeping that search efficient is what separates a toy simulation from a production fluid engine. Master the spatial hash, tune the smoothing length, and you hold one of the most versatile simulation tools in computational science.