dulcimer#

One plucked string and the box it sounds over, simulated as two coupled wave equations and solved live in the browser. Both solvers are written as MATLAB, run by numbl and compiled to WebGPU compute shaders; the pressure field is drawn straight out of the buffer the solver writes, a microphone records the air at one point every timestep, and the recording plays back at the pitch a microphone there would have heard.

This is the instrument-shaped sibling of acoustic-scattering-2d, whose MATLAB-to-WGSL compiler pipeline it reuses, and of acoustic-scattering-3d, whose volume renderer it adapts. What is new here is that two equations run together, on two different grids, with the coupling between them part of the editable model.

The equations#

The string carries transverse displacement u (metres) on a line of nodes:

u_tt = cs^2 u_xx - kap^2 u_xxxx - 2 sig0 u_t + 2 sig1 (u_xx)_t

cs = 2 Ls f0 is the wave speed a string of length Ls needs to sound the fundamental f0. The fourth-derivative term is bending stiffness, expressed through the inharmonicity coefficient B (partial n lands near n f0 sqrt(1 + B n^2)); sig0 is plain decay, quoted as a 60 dB time; sig1 damps high frequencies faster than low, which is why a plucked note starts bright and mellows. The pluck is an initial condition: a triangle drawn to the pluck point, released from rest. This is the standard stiff-string formulation (Bilbao, Numerical Sound Synthesis, ch. 7), with the standard explicit scheme.

The air carries pressure p on a rectangular grid around the instrument:

p_tt + 2 sig p_t = c^2 lapw(p, wall) + s

lapw is a wall-masked Laplacian, and it is how the body is rigid. Each face of the 7-point stencil is scaled by the mask at the neighbour it reads, so a face into the shell carries no flux: the discrete Neumann (sound-hard) condition, in the divergence form div(w grad p). The flat siblings make walls out of fast material instead, which is unavailable here: wood at about 4000 m/s would cut the global timestep twelve-fold, where a mask costs nothing. What a mask cannot do is absorb, so the scene paints a thin lossy skin over the shell's surfaces, and that parameter is what tames the cavity's ring.

The coupling from string to air is one-way and comes in two routes, each with its own gain. string radiation injects the string's acceleration along its own line of cells: the direct route, idealized, since a thin string is in reality a very poor radiator. bridge drive takes the string's pull at its bridge end (the tension times the arriving slope) and drives a patch of air just above the top plate: the real instrument's main route, minus the plate's own resonances, which this model does not carry. Either gain at zero switches that route off. The air does not push back on the string; at these amplitudes the back-reaction is far below everything else the model already neglects.

Two files, two grids#

As in the 2d sibling, the model (init and step, both grids together) is compiled: numbl lowers it to typed IR and this project's backend emits one WGSL kernel per source line, with five host-provided operations — dxx, dxxxx on the string, lapw in the air, and spread/bridge carrying the string into the air — as the only places anything reads a neighbour. The scene (the body: walls, absorption, and the two coupling profiles) is interpreted on the CPU once per edit, so it has the whole MATLAB subset available and costs no recompile.

The two grids are not independent. The air's CFL condition fixes the timestep, and the string then chooses the finest node spacing that is stable at that dt for anything the sliders can reach, so moving a parameter never forces a recompile. In numbers, at the default 128-cell grid: 7.8 mm cells, dt = 6.6 µs (a 152 kHz sample rate), and a string of 59 nodes.

Listening#

The microphone is a one-thread GPU dispatch riding in the same submission as each timestep, so recording costs nothing and no readback happens until you ask to listen. Because everything is in SI units, playback needs no translation: one sample per timestep at 1/dt per second is real time at real pitch. A .wav download (48 kHz, resampled) is one button over.

Two ways to run. Watching takes a few timesteps per display frame, the wave crawling in slow motion. Render note runs the solver flat out with no display until the requested seconds of audio exist, then plays them. On a discrete GPU a second of audio takes very roughly half a minute to render at the default grid; on an integrated one it can take several minutes, which is what the draft 64 grid is for (honest only to about 2.7 kHz, but several times faster).

Honest limitations#

The grid resolves sound to about 5.5 kHz at the default size, so the top of the timbre is simply absent, and the note is duller than a real instrument. The soundboard has no modes of its own: the bridge drives the air directly, so the body colours the sound only through its cavity and hole. The wall mask is perfectly rigid apart from its absorption skin, and the box walls are drawn at grid resolution, so plates thinner than about two cells leak. The coupling gains are physically arbitrary; the equations are linear, so they set relative balance, not absolute loudness.

Running it#

npm install
npm run dev       # local dev server
npm test          # solver checks against desktop WebGPU (Google Dawn)
npm run smoke     # headless-browser check of the built page (environment permitting)

npm test compiles the actual .m files to actual shaders and checks physics: the measured pitch of the plucked string, the 60 dB decay time, that a sealed box keeps sound out and a sound hole lets the cavity speak, and that the microphone's spectrum sits on the string's partial comb.