Differential Equations in Simulations — ODEs, PDEs and Numerical Methods Explained

What is a Differential Equation?

A differential equation is an equation that relates a function to its own derivatives. In physics, the "function" is usually position, velocity, temperature or chemical concentration, and the "derivative" is its rate of change over time or space.

An ODE (ordinary differential equation) involves derivatives with respect to a single variable, usually time. A PDE (partial differential equation) involves derivatives with respect to multiple variables — time and space simultaneously.

Part 1 — Ordinary Differential Equations (ODEs)

The Simplest ODE: Exponential Growth and Decay

The ODE dy/dt = ky has the exact solution y(t) = y₀ · e^kt. When k < 0 it describes radioactive decay, drug clearance and population decline. When k > 0 it describes compound interest and bacterial growth. Every numerical method is just an approximation to this exact solution.

Euler Method — Fast but Inaccurate

The Euler method advances the state by one small time step Δt using only the current derivative (slope). It is the simplest possible integrator:

          // Euler method — O(Δt) global error
          function
          eulerStep(state, dt) {
          const deriv =
          derivatives(state);
          // e.g. { vel: acc, pos: vel }
          return { pos: state.pos + deriv.vel * dt, vel:
          state.vel + deriv.acc * dt }; }
        

The problem: error accumulates as O(Δt) per step. For oscillating systems like a pendulum, Euler slowly adds energy each step — the pendulum swings wider and wider until it escapes. You can see this in the Pendulum simulation by enabling "Euler" mode.

Runge-Kutta 4 (RK4) — The Workhorse Integrator

RK4 evaluates derivatives at four points within the time step and combines them with a weighted average. Error drops to O(Δt⁴) — 10,000× more accurate than Euler for the same step size.

          // RK4 — O(Δt⁴) global error, ~4× more expensive than Euler
          function
          rk4Step(state, dt) {
          const k1 =
          derivatives(state);
          const k2 =
          derivatives(add(state,
          scale(k1, dt/2)));
          const k3 =
          derivatives(add(state,
          scale(k2, dt/2)));
          const k4 =
          derivatives(add(state,
          scale(k3, dt)));

          // Weighted average: 1/6 k1 + 1/3 k2 + 1/3 k3 + 1/6 k4
          return add(state,
          scale( addAll(k1,
          scale(k2, 2),
          scale(k3, 2), k4), dt
          / 6
          )); }
        

Verlet Integration — Energy-Conserving for Mechanics

Verlet integration is the preferred method for particle simulations. It is symplectic — it conserves a modified energy over long integrations — which is why cloth, soft-body and molecular dynamics simulations use it instead of RK4:

          // Velocity Verlet — symplectic, O(Δt²) but energy-stable
          pos_new = pos + vel * dt + 0.5 * acc * dt *
          dt; vel_new = vel + 0.5 * (acc + acc_new) *
          dt;
        

Part 2 — Partial Differential Equations (PDEs)

The Heat Equation

The heat equation describes how temperature diffuses through a medium over time: ∂T/∂t = α ∇²T. The spatial Laplacian ∇²T (second derivative in space) is discretised with finite differences on a grid.

2D heat equation — finite difference explicit scheme

T[i][j] += α · dt / dx² · (T[i+1][j] + T[i-1][j] + T[i][j+1] + T[i][j-1] − 4·T[i][j])

Stability condition: α · dt / dx² ≤ 0.25

The Wave Equation

Sound, water waves and strings all obey ∂²u/∂t² = c² ∇²u. Unlike the heat equation (first time derivative — "memory-less"), the wave equation's second time derivative means disturbances travel rather than diffuse. The Chladni figures simulation solves the 2D wave equation on a square plate.

          // 2D wave equation — explicit finite difference (bath waves,
            Chladni)
          for (let i =
          1; i < N-1; i++)
          { for (let j =
          1; j < N-1; j++)
          { const laplacian = u_prev[i+1][j] + u_prev[i-1][j] + u_prev[i][j+1] + u_prev[i][j-1] -
          4 * u_prev[i][j]; u_next[i][j] =
          2 * u_curr[i][j] - u_old[i][j] + c2 * dt2 *
          laplacian; } }
        

Reaction-Diffusion PDEs

The Gray-Scott system (seen in the Chemistry Spotlight) is a PDE system where diffusion (Laplacian terms) combines with nonlinear reaction terms. The spatial Laplacian is computed the same way as the heat equation — it is the reaction term that creates the Turing instability and pattern formation.

Choosing the Right Method

Method	Error order	Cost	Best for
Euler	O(Δt)	1 eval/step	Prototyping; never production physics
Verlet	O(Δt²)	1–2 evals/step	Particle systems, molecular dynamics, cloth
RK4	O(Δt⁴)	4 evals/step	Orbital mechanics, ODEs needing accuracy
FDM (explicit)	O(Δt, Δx²)	O(N²) per step	Heat eq., wave eq. on fixed grids
FDM (implicit)	O(Δt², Δx²)	Solve linear system	Stiff PDEs, large time steps needed

Where These Methods Appear on the Site

Binary Stars / Solar Sail

RK4 integrates gravitational ODEs. Energy drift over 10,000 orbits is <0.01% — impossible with Euler.

Cloth Simulation

Verlet integration over a spring-mass grid. Position-based constraints handle collisions without velocity corrections.

Bath Waves

Explicit FDM wave equation on a 256×256 grid. Stability enforced by the CFL condition: c·Δt/Δx ≤ 1.

Reaction-Diffusion

Gray-Scott coupled PDEs on a 512×512 texture. Finite differences run in a GLSL fragment shader for GPU parallelism.

Drug Diffusion

Fick's second law (heat equation variant) solved with explicit FDM. Stability condition: D·Δt/Δx² ≤ 0.5.

Chladni Figures

2D wave equation eigenvalue problem solved analytically for mode shapes, then rendered with particle deposition on nodal lines.

The key insight: Every numerical method is an approximation to the same underlying differential equation. The question is always: how much accuracy do you need, and how much computation can you afford? RK4 is not always the answer — for long particle simulations, Verlet's symplectic property outperforms RK4's accuracy.