Barnes-Hut N-Body

Imagine trying to simulate a galaxy. Every star exerts a gravitational pull on every other star. With $n$ stars, that is $n \times (n-1)$ pairs — roughly $n^2$ force calculations per time step. At one million stars you need a trillion calculations every frame. At a billion stars the number becomes astronomical in both senses of the word.

This is the N-body problem: given $n$ point masses, compute the net gravitational force on each one so you can advance their positions by a tiny time step and repeat. It underlies everything from planetary mechanics to cosmological simulations.

The naive $O(n^2)$ algorithm is exact but hopeless at scale. In 1986 Josh Barnes and Piet Hut published a clever approximation: instead of summing every pairwise interaction, group distant particles together and treat the group as a single mass at its center of mass. Their algorithm runs in $O(n \log n)$ — fast enough to simulate millions of bodies on ordinary hardware.

The key data structure is a quadtree (in 2-D) or octree (in 3-D): a recursive spatial subdivision that lets you decide, for each particle, which distant clusters are far enough away to be summarized safely.

The simulation below spawns particles in a rotating disk and evolves them under gravity. Toggle Show quadtree to see how the algorithm partitions space — each cell summarizes all its particles as a single mass when a particle is far enough away.

<div class="controls">
  <label>{{lbl_theta}} <input id="theta" type="range" min="0.1" max="1.5" step="0.05" value="0.5"> <span id="theta-val">0.50</span></label>
  <label><input id="show-tree" type="checkbox"> {{lbl_show_tree}}</label>
  <button id="btn-pause" type="button">{{btn_pause}}</button>
  <button id="btn-reset" type="button" class="ghost">{{btn_reset}}</button>
</div>
<canvas id="canvas" width="480" height="320"></canvas>
<div class="stats" id="stats"></div>

* { box-sizing: border-box; }
body { margin: 0; background: #0b0d14; font-family: system-ui, sans-serif; color: #c8d4e0; }
.controls { display: flex; align-items: center; gap: .8rem; flex-wrap: wrap; padding: .5rem .6rem; background: #131622; border-bottom: 1px solid #1e2538; font-size: .85rem; }
label { display: flex; align-items: center; gap: .3rem; cursor: pointer; user-select: none; }
input[type=range] { width: 90px; accent-color: #5b9bd5; }
button { font: 600 13px system-ui; padding: .3rem .8rem; border-radius: 6px; cursor: pointer; border: 1px solid #5b9bd5; background: #5b9bd5; color: #fff; }
button.ghost { background: transparent; color: #5b9bd5; }
canvas { display: block; width: 100%; max-width: 480px; }
.stats { font-size: .78rem; padding: .35rem .6rem; background: #131622; color: #7a91a8; min-height: 1.6em; }

// Code not found

Notice how the quadtree cells near a particle are small and numerous (high resolution), while distant regions are large single boxes (one interaction instead of hundreds). The theta parameter $\theta$ controls the trade-off: smaller values are more accurate but slower; larger values are faster but rougher. The default $\theta = 0.5$ is the value Barnes and Hut recommended in their original paper.

Naive N-body is $O(n^2)$ per time step: every particle loops over every other particle. Doubling $n$ quadruples the work.

Barnes-Hut breaks that with two observations:

1. Build a quadtree. Recursively subdivide the simulation region into four quadrants until each leaf holds at most one particle. This costs $O(n \log n)$ in the average case.
2. Apply the opening-angle criterion. For each particle $p$, walk the tree. At each internal node (a cell of width $s$ at distance $d$ from $p$), check whether $s / d < \theta$. If yes, the whole cell's mass can stand in for all its particles — one interaction instead of potentially thousands. If no, recurse into the children.

The fraction of nodes that survive the opening-angle test is $O(\log n)$ on average, so the total cost per particle is $O(\log n)$, giving the full step $O(n \log n)$.

The approximation error in the force on particle $p$ is bounded by a term proportional to $\theta^2 \cdot (s/d)^2$, so halving $\theta$ roughly quarters the error while increasing runtime. In practice $\theta \in [0.5, 1.0]$ gives results close enough for astrophysical simulations.

Compare this to the fast Fourier transform: both algorithms replace a seemingly irreducible $O(n^2)$ task with $O(n \log n)$ by exploiting structure in the problem — spatial proximity here, frequency decomposition there.

Anywhere you sum a long-range interaction over many particles, Barnes-Hut or one of its cousins shows up:

• Astrophysical simulations: galaxy formation, dark-matter halos, stellar clusters — Barnes-Hut (and its successor the fast multipole method) made large-scale cosmological simulations feasible on real machines.
• Molecular dynamics: electrostatic forces between atoms follow the same $1/r^2$ law as gravity. The particle-mesh Ewald and fast multipole methods are direct descendants of the same idea.
• Force-directed graph layouts: tools like D3.js use Barnes-Hut to repel graph nodes, turning a layout that would freeze on large graphs into one that runs in the browser.
• Smoothed-particle hydrodynamics (SPH): fluid simulation where each "particle" carries a small blob of fluid; neighbor queries use the same spatial tree.
• Machine learning: $t$-SNE uses Barnes-Hut to approximate gradient repulsions among $n$ data points, enabling visualization of datasets with tens of thousands of items.

The common thread is long-range interaction + spatial locality: nearby things need full precision, distant things can be grouped. Whenever those two conditions hold, a spatial tree can replace $O(n^2)$ with $O(n \log n)$.

To see a related $O(n \log n)$ idea in a completely different domain, compare with the fast Fourier transform.

Barnes and Hut's insight was not to compute less — it was to compute smartly. A distant cluster of a thousand stars can be replaced by one number (its total mass) and one point (its center of mass) without meaningfully affecting the trajectory of a star far away. The quadtree automates the decision of what counts as "far away."

The result is an algorithm that scales to millions of particles on a laptop, to billions on a supercomputer — all because of a single ratio check: $s / d < \theta$.

That ratio check is an example of a broader pattern in algorithm design: controlled approximation. Rather than demanding exact answers at every step, accept a bounded error and win orders of magnitude in speed. The same philosophy drives fast Fourier transform algorithms, compressed sensing, and the approximate nearest-neighbor searches that power modern recommendation systems.

Gravity is, in the end, a tractable problem — not because the physics simplifies, but because the geometry does.