Tip: Visualising ODEs and PDEs in the Browser

Differential equations describe how things change — a planet's velocity, the temperature across a metal plate, the swirl of a fluid. The browser cannot solve them symbolically in real time, but it does not need to. It just needs to advance the state by a small time step, draw the result, and repeat. Get that loop right and a static equation turns into a living, pannable, parameter-tweakable picture. Get it wrong and your simulation either crawls or explodes into NaN within a few seconds.

Below is the workflow we use across the catalogue, from ordinary differential equations (ODEs) that track a handful of variables to partial differential equations (PDEs) that evolve a whole grid of values.

1. Pick the Right Integrator: Euler vs RK4

An ODE solver answers one question repeatedly: given the current state and the rate of change, what is the state a tiny moment later? The simplest answer is the forward Euler method — take the current slope and step straight along it.

// state = {x, v}, derivs(s) returns {dx, dv}
function eulerStep(s, dt) {
  const d = derivs(s);
  return { x: s.x + d.dx * dt, v: s.v + d.dv * dt };
}

Euler is cheap and perfectly adequate for gentle, slow systems or particle clouds where a little drift is invisible. But it accumulates error fast, and for oscillating or chaotic systems it injects fake energy until the orbit spirals out. That is where fourth-order Runge–Kutta (RK4) earns its keep. RK4 samples the slope four times per step — once at the start, twice in the middle, once at the end — and blends them into a far more accurate estimate.

function rk4Step(s, dt) {
  const k1 = derivs(s);
  const k2 = derivs(add(s, scale(k1, dt / 2)));
  const k3 = derivs(add(s, scale(k2, dt / 2)));
  const k4 = derivs(add(s, scale(k3, dt)));
  return add(s, scale(add4(k1, k2, k2, k3, k3, k4), dt / 6));
}

RK4 costs four derivative evaluations instead of one, but it lets you take much larger time steps for the same accuracy — usually a net win. Our Lorenz attractor in 3D relies on RK4: the butterfly shape only stays crisp because the integrator does not bleed energy. For stiff or energy-conserving mechanical systems, a symplectic variant like velocity Verlet is often the better pick, but RK4 is the safe default when you are not sure.

2. Finite-Difference Grids for PDEs

PDEs evolve a field — a value at every point in space — rather than a few scalars. The standard browser-friendly approach is the finite-difference method: store the field as a flat Float32Array indexed as i = y * width + x, and approximate spatial derivatives by comparing each cell to its neighbours.

The workhorse is the discrete Laplacian, which drives diffusion, heat flow, and wave propagation:

// 5-point Laplacian on a width×height grid
const i = y * width + x;
const lap = field[i - 1] + field[i + 1]
          + field[i - width] + field[i + width]
          - 4 * field[i];
next[i] = field[i] + alpha * lap * dt;

Always write into a separate output buffer and swap the two arrays after the full sweep — updating in place corrupts the neighbours you have not visited yet. This double-buffer pattern is exactly what powers grid simulations like plasma instability, where a charged fluid field is marched forward cell by cell, and the diffusion-style smoothing behind granular heating, where local collisions spread energy through a packed bed of grains.

3. Respect the Stability Condition

Here is the trap that catches everyone the first time: explicit finite-difference schemes are only stable if the time step is small enough relative to the grid spacing. For a diffusion equation in two dimensions the CFL (Courant–Friedrichs–Lewy) condition roughly requires alpha * dt / dx² < 0.25. Push past it and the field does not just get inaccurate — it oscillates wildly and overflows to infinity within a handful of frames.

Two practical fixes keep you safe without slowing the visible animation. First, run several small physics sub-steps per rendered frame — the maths stays stable while the screen still updates 60 times a second. Second, clamp or gently damp the field each step so a stray spike cannot snowball. Wave-style PDEs, like the dispersive ripple pattern in our Kelvin wake simulation, are especially sensitive here: the V-shaped wake only holds its shape when the wave speed and time step stay in balance.

4. Decouple Physics from Rendering

A common mistake is tying the integration step directly to requestAnimationFrame and passing the real elapsed time as dt. When the tab stalls, that delta balloons to hundreds of milliseconds and a single giant step detonates your stable scheme. The fix is a fixed-timestep accumulator: advance the physics in constant slices and only render once per frame.

let acc = 0;
const FIXED = 1 / 120; // physics seconds per sub-step

function frame(now) {
  acc += Math.min((now - last) / 1000, 0.1); // cap the spike
  last = now;
  while (acc >= FIXED) { step(FIXED); acc -= FIXED; }
  render();
  requestAnimationFrame(frame);
}

Capping the accumulated time (Math.min(..., 0.1)) prevents the dreaded "spiral of death" where a slow frame demands more sub-steps than the next frame can afford. This single pattern is the difference between a simulation that survives a background tab and one that returns to a screen full of garbage.

5. Render Fast: Batch Your Canvas Calls

Even perfect maths feels broken at 20 FPS. Canvas 2D is plenty fast for these visualisations if you stop fighting it:

• Draw fields with putImageData. For a full grid, write RGBA bytes straight into an ImageData buffer and blit the whole frame in one call — never one fillRect per cell.
• Batch particle paths. For ODE swarms, open a single beginPath(), push every moveTo/lineTo, and issue one stroke(). Each separate stroke is a costly state change.
• Avoid per-frame allocation. Reuse your typed arrays and state objects. Garbage-collection pauses are the most common cause of visible 60 FPS stutter.
• Match the device pixel ratio once. Scale the canvas backing store to devicePixelRatio at setup, not every frame, for crisp lines on retina displays without per-frame overhead.

Together these turn a sluggish prototype into something that feels instantaneous, even on a mid-range phone.

Putting It Together

The recipe is the same whether you are evolving three coupled ODEs or a million-cell PDE field: choose an integrator that matches how sensitive your system is (Euler for the forgiving, RK4 for the rest), lay your field out as a flat typed array with a double-buffered finite-difference update, keep your time step under the stability limit, decouple physics from rendering with a fixed accumulator, and batch every Canvas call you can. Master those five habits and almost any differential equation in a textbook becomes a simulation you can ship.

Want to see these techniques in motion? Explore the Lorenz attractor, plasma instability, Kelvin wake, and granular heating simulations — every one of them is a differential equation running the loop described above. Or browse the full simulation catalogue to find the equation behind your favourite phenomenon.

Related Posts