N-body Simulation
100,000+ gravitating bodies at interactive frame rates, entirely in the browser.
Direct-sum gravity is O(N²). In JavaScript that caps out in the low thousands of bodies - far short of anything that looks like a galaxy.
Three interchangeable engines: brute-force CPU, a Barnes–Hut quadtree, and a WebGPU compute-shader kernel. Ping-pong storage buffers let the vertex shader read simulation output directly instead of round-tripping to CPU RAM, the CPU engines use SharedArrayBuffer for zero-copy worker IPC, and a symplectic leapfrog integrator keeps orbits stable across thousands of periods.
Over 100,000 bodies at interactive frame rates in a browser tab, with an isothermal dark-matter halo and Salpeter-IMF stellar masses. No server, no install.
What is this?
An N-body simulation works out how a group of objects move when each one pulls on all the others through gravity. It's the standard tool astronomers use to study how galaxies form and change.
What makes it hard is the amount of computation. Every star pulls on every other star, so 10,000 stars means about 100 million pair calculations for a single step forward in time, and a smooth animation needs dozens of steps per second. This project makes that work in an ordinary web browser by offering four "engines" that solve the same problem in different ways, from a simple exact method up to one that runs on the graphics card.
At a glance
- Fast: 200,000 stars at 60 frames per second on a laptop's built-in graphics chip.
- Checked against real physics: 22 automated test files compare the simulation with known exact answers, such as an orbit that must return to its starting point. They run on every code change, and the site only deploys if they pass.
- Reproducible: every simulation starts from a random seed, so a shared link rebuilds exactly the same galaxy.
- Measured claims: a built-in benchmark records the speed numbers below, and a live panel shows how well the simulation conserves energy.
Simulation Presets
- Galaxy (spiral). A flat, spinning disk of stars inside an invisible dark-matter halo. The stars' own gravity holds the disk together, and spiral arms keep forming, fading and re-forming. The disk is set up to be only just stable: stable enough not to collapse into clumps, loose enough to grow spiral arms. Astronomers measure this balance with the Toomre Q number, which is set to 1.3 here.
- Accretion disk. A supermassive black hole at the centre with a disk of lighter stars orbiting it, like planets around the Sun. The black hole dominates, so the stars follow clean elliptical orbits (Keplerian orbits). Gas friction is not modelled, so nothing spirals inwards. The timestep shrinks automatically near the black hole, where orbits are fastest. Dark matter is off by default.
Features
Four physics engines
All four calculate the same gravity and can be switched while the simulation is running. Each one has a star limit based on what it can sustain smoothly, and the interface tells you when you hit it.
| Engine | How it works | Where it runs | Max stars |
|---|---|---|---|
| Brute Force | Exact: every star against every other star, $O(N^2)$ | Main thread | 20 000 |
| Barnes-Hut | Approximate: groups distant stars together, $O(N \log N)$ | Main thread | 50 000 |
| Worker | Barnes-Hut on a background thread, so the page stays responsive | Web Worker, sharing memory with the page | 50 000 |
| WebGPU | Exact, with thousands of stars calculated in parallel | Graphics card (GPU) | 200 000 |
WebGPU is the default. If the browser has no usable GPU, or the GPU stops responding mid-run, the simulation switches to Barnes-Hut. The Worker engine needs a page setting called cross-origin isolation and falls back the same way without it.
Physics
- Leapfrog integration. This is how the simulation steps forward in time. Unlike simpler methods, it keeps the total energy stable over thousands of orbits instead of letting it creep up or down. Positions and velocities are updated half a step apart, which gives it that stability.
- A galaxy that starts in balance. Each star's starting speed is calculated from the gravity the actual generated stars produce, not from an idealised formula. Starting speeds are then randomised just enough to hit the Toomre Q target. The result is a disk that holds its shape from the first frame instead of collapsing or flying apart.
- Adaptive timestep. The step size shrinks when stars orbit fast or pass close to each other, so fast motion isn't skipped over. A minimum size stops the simulation from grinding to a halt.
- Dark matter halo. It's modelled as a smooth background pull rather than millions of extra particles. This is enough to give the galaxy its flat rotation curve, meaning outer stars orbit about as fast as inner ones, as in real galaxies.
- Realistic star masses and colours. Star masses follow the Salpeter distribution (many small stars, a few heavy ones), and colours follow the main sequence of the Hertzsprung-Russell diagram (heavy stars blue-white, light stars orange-red).
- Active and passive stars. To improve performance, only a subset of stars (the active ones) produce gravity, but every star feels it. The active set is always "the first
activeCountstars", and every engine applies that same rule, so all four simulate exactly the same system. In the accretion preset the active stars are the heaviest ones. In the galaxy preset they are a random sample that carries the disk's full mass.
Diagnostics
- Live energy panel. In a correct simulation total energy should stay almost constant, so its drift is the best single sign of accuracy. The panel plots that drift ($\Delta E/E_0$) for the active stars. The exact energy calculation is expensive, so it's spread over several frames to avoid stutter. The panel is hidden on the WebGPU engine, because the star data never leaves the GPU there.
- Share links. The Share button copies a link that rebuilds the same starting galaxy. See Reproducibility.
- Built-in benchmark. Add
?benchto the URL to time every engine at every star count and get a ready-made results table. - Accurate speed readout. The top bar shows pair calculations per second, counted exactly for each engine so engines can be compared fairly. The GPU time comes from the GPU's own timer where the hardware supports it.
Key Decisions
| Decision | Why |
|---|---|
| WebGPU instead of WebGL | WebGPU has proper compute shaders and storage buffers. WebGL was built for drawing, not general calculation, so the parallel gravity calculation fits WebGPU much better. |
| Tiled GPU kernel | The default GPU kernel loads stars in batches of 64 into fast on-chip memory that a group of threads shares, instead of each thread fetching every star from slower main GPU memory. The simple (naive) kernel still ships alongside it, and the benchmark checks that both give the same results, so the speed-up is verified rather than assumed. |
| Shared memory for the Worker engine | The worker runs Barnes-Hut on a background thread. Sharing one block of memory (SharedArrayBuffer) with the page avoids copying every star's position between threads on every frame. Only the Worker engine uses it. |
| One force function for all CPU engines | Brute force calculates each pair of stars twice instead of using Newton's third law to halve the work. That costs 2× the arithmetic, but every CPU engine then uses the same tested force function (pairwiseAccel). As a result, any difference between Barnes-Hut and brute force can only come from Barnes-Hut's approximation. |
| Active stars chosen by position in the list | Every engine treats stars [start, activeCount) as the gravity sources, set once by the preset. Choosing them by mass instead would let engines disagree. That bug really happened once, and it's described in Challenges. |
| Two GPU buffers, swapped each step | The shader reads from one buffer and writes to the other, so no thread ever reads a value another thread is overwriting. The drawing code reads the output buffer directly on the GPU, with no copy back to the CPU. |
| Leapfrog integrator | Picked over Euler and Runge-Kutta because it keeps energy stable over thousands of orbits, where those methods slowly drift. |
| Energy panel measures only the active stars | Passive stars feel gravity but don't produce it. That one-way pull saves a lot of work, but it means the energy of the whole system is not expected to stay constant, so measuring it would only show a side effect of the shortcut. The active stars on their own form a proper closed system, so their energy is the meaningful thing to track. |
Challenges & Lessons
Keeping two threads in step
The Worker engine coordinates the page and the background thread with Atomics.wait and Atomics.notify on a single shared status flag.
- The hard part wasn't the handshake, it was the bookkeeping. If the worker is still busy when a new frame arrives, that frame's step is skipped. The simulation clock may only advance for steps that actually finished, or it quietly runs too fast.
- The order of operations matters. Drawing waits on the same flag, so the page never draws half-updated positions. When I started the next step before drawing the last one, the canvas froze, while every counter still reported normal progress.
Garbage collection slowing every frame
Barnes-Hut rebuilds its tree every frame, and creating thousands of new objects each time made the browser's garbage collector pause constantly. I fixed it in two ways:
- The tree reuses its nodes from a pool instead of creating new ones.
- The force functions write into an object the caller passes in, instead of returning a new object on every call.
Learning to think like a GPU
On a GPU, groups of threads run in lockstep, so an if that sends threads different ways makes them wait for each other. The tiled kernel also has sync points (workgroupBarrier) that every thread in the group must reach. As a result, a thread with no star to process can't simply exit early. It carries a "valid" flag and runs to the end with everyone else. The classic bug in tiled N-body code is a barrier that some threads never reach.
Four engines that quietly disagreed
- The bug. Brute force and the GPU chose gravity sources by their position in the list. Barnes-Hut chose them by mass, keeping stars above a threshold.
- Why it hid. In the accretion preset, stars are sorted by mass, so both rules picked the same stars. In the galaxy preset every disk star has the same mass, so Barnes-Hut counted far more of them as sources. It simulated a disk 3.3× heavier than the one the starting conditions were balanced for, and its energy drift was 100,000 times worse than brute force.
- Why the tests missed it. Every galaxy test used so few stars that all of them were active, so the two rules never differed.
- How it was found and fixed. A test comparing Barnes-Hut (with approximation off) against brute force, with the active/passive split switched on, exposed it in one step. That comparison is now a permanent test, and every engine reads the same
activeCountvalue. - The lesson. When several engines must agree on a rule, give them one shared value to read rather than letting each one re-implement the rule.