The Starting Point: What Is an EM Wave?
James Clerk Maxwell unified electricity and magnetism in 1865 with four equations. The key insight for this simulation is the wave equation that falls out of them: a changing electric field induces a magnetic field, and vice versa. The two fields chase each other through space at the speed of light.
In 2D, we work with the transverse magnetic (TM) mode: the electric
field has only a z-component Ez, and the magnetic field
has x and y components Hx and Hy. This
reduces Maxwell's curl equations to three update equations — simple
enough to evaluate on every grid cell every frame.
FDTD: The Yee Lattice
The Finite Difference Time Domain (FDTD) method, introduced by Kane Yee in 1966, discretises Maxwell's equations on a staggered space-time grid. Electric and magnetic field components are offset by half a spatial and temporal step, which automatically satisfies Gauss's law (∇·B = 0) without any correction step.
// TM-mode FDTD update equations (2D, uniform grid)
// E and H are interleaved: H updates first, then E
// Update Hx and Hy at time step n+½
Hx[i][j] -= (Cz / dy) * (Ez[i][j+1] - Ez[i][j])
Hy[i][j] += (Cz / dx) * (Ez[i+1][j] - Ez[i][j])
// Update Ez at time step n+1
Ez[i][j] += Cz * (
(Hy[i][j] - Hy[i-1][j]) / dx -
(Hx[i][j] - Hx[i][j-1]) / dy
)
// Cz = dt / (ε₀ · μ₀) — Courant number
The Courant stability condition requires
dt ≤ dx / (c · √2) for a 2D uniform grid. We set
dt = 0.9 · dx / (c · √2) — 90 % of the stability limit
for a small safety margin. With a 200×200 grid and dx = 0.01 m, the
time step is about 21 ps per update.
The Build Timeline
Started with a 100×100 JavaScript array, pure CPU, single point source. The wave propagated correctly — circular rings spreading outward. But it reflected off the edges perfectly, so the whole domain filled up with standing waves within seconds. Clearly needed absorbing boundaries.
PML is the standard FDTD absorbing boundary: a lossy sponge region around the domain where the fields are attenuated without reflection. Implementing it doubled the code complexity but eliminated reflections almost perfectly. The domain now felt "open" — waves exit cleanly.
The CPU implementation ran at ~30 FPS on a 100×100 grid. Extending to 200×200 (the minimum grid for useful visualisation) dropped it to 8 FPS. The FDTD update equations are embarrassingly parallel — every cell can be updated independently. Moved the update loop to a WebGL fragment shader writing into a floating-point texture pair (ping-pong scheme). Jumped to 60 FPS at 400×400.
Mapped the Ez field to a blue-red colour scale (negative to positive), added a sinusoidal point source modulated by a Gaussian envelope to avoid the initial "click", and implemented a double-slit aperture for the classroom-famous diffraction pattern. Added dipole source, plane wave source, and a movable dielectric region the user can place on the canvas.
Overlaying Hx/Hy as arrows using instanced line segments in WebGL. Arrows are subsampled to a 20×20 grid, scaled by field magnitude, and rotated to field direction. The combination of Ez colour map and H-field arrows gives a visceral sense of how the two fields are in quadrature.
The Hardest Part: Floating-Point Textures
WebGL1 does not support reading from a floating-point texture that
you've just written to — the "render to float texture" extension
(OES_texture_float) exists but isn't universally
writable. The solution is a ping-pong buffer: two
framebuffers, swapping roles each frame. The shader reads from buffer
A and writes to buffer B; next frame it reads from B and writes to A.
On iOS Safari the extension wasn't available at all, requiring a fallback to 16-bit half-float textures with reduced dynamic range. Field values are normalised before storage and rescaled on read — a small precision loss, but invisible at the colour resolution of the display.
The Dielectric Permittivity Texture
To let users place materials (glass, metal, etc.) in the domain, the FDTD update coefficients need to vary spatially. We bake the permittivity εr into a separate R8 texture and sample it during the update shader. Pressing a toolbar button sets a "brush mode" — subsequent mouse events paint permittivity values into the material texture via a small secondary framebuffer.
// In GLSL update shader — reading material texture float eps_r = texture2D(u_material, vTexCoord).r * 10.0 + 1.0; // Range: [1.0, 11.0] — air to glass-like dielectric float Cz_local = dt / (eps_r * eps0 * mu0 * dx * dx); // Same update formula, but Cz varies with material
Classroom Use: What Works and What Doesn't
We got feedback from several physics teachers after the simulation launched. The features that proved most valuable in the classroom:
- Double-slit diffraction preset — students can directly observe Huygens' principle: each slit becomes a secondary source, and the interference pattern emerges from the FDTD physics without any hand-tuning.
- Speed of light slider — scaling c down by 10× makes wave propagation slow enough to trace with the eye. Making it visceral.
- Polarisation mode toggle — switching between TE and TM modes shows how the two polarisations of light are just rotated versions of the same field structure.
What didn't work: the initial free-draw dielectric tool confused students because "drawing glass" looks just like "drawing air" before the wave hits. We added a persistent colour overlay for placed materials (amber tint = dielectric, grey = conductor) which helped enormously.
Try it yourself: Open the EM Waves simulation, hit the "Double Slit" preset, and drag the slit separation slider. You'll see the diffraction pattern change in real time — the same physics as Thomas Young's 1801 experiment, running in your browser.
Next Steps for the Simulation
The 2D TM-mode FDTD is a strong teaching tool, but there are natural extensions I want to build:
- 3D mode — all six field components on a 3D Yee lattice. Requires WebGPU compute for the update pass; the data doesn't fit in a 2D texture at useful grid sizes.
- Dispersive materials — Drude model for metals (frequency-dependent permittivity), enabling plasmonic resonance visualisation.
- Near-to-far-field transform — compute the far-field radiation pattern (antenna pattern) from the near-field stored on a Huygens box around the source.