d8098edAn interactive tutorial on the wave equation in acoustics: strings, drums and rooms, measured by earJeremy Magland 1import type { Chapter } from '../../ui/types.ts'
2import type { Num, Params } from '../../engine/types.ts'
3import driver from './driver.m?raw'
4import script from './leapfrog.m?raw'
5import text from './text.md?raw'
6import Listen, { type Clip } from '../../audio/Listen.tsx'
7import Field from '../../plot/Field.tsx'
8import { Row, SeriesPlot, fmtCents, fmtErr, fmtNum } from '../common/panels.tsx'
9import { familyColor, familyShapes, ink, method } from '../../plot/palette.ts'
10import { FS } from '../shared.ts'
12interface RunOut {
13 ms: number
14 lam: number
15 N: number
16 unstable: boolean
17 blowup: Num
18 tp: number[]
19 err: Num[]
20 relerr: number
21 tw: number[]
22 yw: Num[]
23 yew: number[]
24 pred1: number
25 pred10: number
26 clipExact: Clip
27 clipFd: Clip
28 clipDiff: Clip
29}
31interface CentsOut {
32 N: number
33 f0: number
34 lams: number[]
35 sel: number[]
36 nall: number[]
37 meas: Num[][]
38 pred: Num[][]
39 ms: number[]
40 agree: number
41}
43interface ClickOut {
44 ms: number
45 N: number
46 lam: number
47 fcut: number
48 rel: number
49 tw1: number[]
50 y1: number[]
51 ye1: number[]
52 tw20: number[]
53 y20: number[]
54 ye20: number[]
55 spec: Num[][]
56 specF: [number, number]
57 specT: [number, number]
58 clipExact: Clip
59 clipFd: Clip
60}
62const F0 = FS / 400 // 110.25 Hz: the pitch at which N = 200 is exactly lam = 1
64function RunPanels({ d }: { d: RunOut }) {
65 return (
66 <>
67 <SeriesPlot
68 title="three periods of the recorded signal, starting one second after the pluck"
69 height={190}
70 xLabel="ms"
71 series={[
72 { x: d.tw, y: d.yew, color: ink.reference, label: "exact (d'Alembert)", width: 2.6 },
73 { x: d.tw, y: d.yw, color: method.fd2, label: `leapfrog, N = ${d.N}, λ = ${fmtNum(d.lam, 3)}`, width: 1.4 },
74 ]}
75 />
76 <SeriesPlot
77 title="relative rms error of the recorded signal, period by period"
78 height={170}
79 yLog
80 yDomain={[1e-16, 10]}
81 xLabel="ms"
82 yLabel="error"
83 series={[{ x: d.tp, y: d.err.map((v) => (v != null && v > 1e-16 ? v : 1e-16)), color: method.fd2, width: 1.6 }]}
84 legend={false}
85 />
86 <Listen
87 items={[
88 { label: 'exact', clip: d.clipExact },
89 { label: 'leapfrog', clip: d.clipFd },
90 { label: 'the difference', clip: d.clipDiff },
91 ]}
92 note={d.unstable ? 'the unstable run is silenced from the moment it blows up' : 'all three at the same scale'}
93 />
94 </>
95 )
96}
98function CentsPanels({ d }: { d: CentsOut }) {
99 const n = d.lams.length
100 const col = (k: number) => familyColor(k, n)
101 const ppwAll = d.nall.map((m) => (2 * d.N) / m)
102 const ppwSel = d.sel.map((m) => (2 * d.N) / m)
103 const measured = (k: number, abs: boolean) => ({
104 y: d.meas[k].map((v) => (v == null ? null : abs ? Math.abs(v) : v)),
105 color: col(k),
106 shape: familyShapes[k],
107 line: false,
108 label: `λ = ${d.lams[k]}`,
109 tip: (i: number) => [
110 `partial ${d.sel[i]}, ${(d.sel[i] * d.f0).toFixed(0)} Hz`,
111 `${fmtNum(ppwSel[i], 3)} points per wavelength`,
112 `measured ${fmtCents(d.meas[k][i])}`,
113 ],
114 })
115 return (
116 <>
117 <SeriesPlot
118 title={`pitch error of every partial on a grid of N = ${d.N}: measured (markers) and the dispersion relation (lines)`}
119 height={260}
120 xLabel="partial frequency, Hz"
121 yLabel="cents"
122 series={d.lams.flatMap((_, k) => [
123 { x: d.nall.map((m) => m * d.f0), y: d.pred[k], color: col(k), width: 1.2, noLegend: true },
124 { x: d.sel.map((m) => m * d.f0), ...measured(k, false) },
125 ])}
126 />
127 <SeriesPlot
128 title="the same errors against points per wavelength: every grid falls on the same curves"
129 height={240}
130 xLog
131 yLog
132 xDomain={[1.8, 500]}
133 yDomain={[1e-4, 1500]}
134 xLabel="points per wavelength"
135 yLabel="|cents|"
136 hlines={[{ y: 5, label: 'a commonly quoted discrimination threshold, 5 cents' }]}
137 refs={[{ p: -2, through: [40, 3], label: 'slope −2', span: 0.4 }]}
138 series={d.lams.flatMap((l, k) =>
139 l === 1
140 ? []
141 : [
142 { x: ppwAll, y: d.pred[k].map((v) => (v == null ? null : Math.abs(v))), color: col(k), width: 1.2, noLegend: true },
143 { x: ppwSel, ...measured(k, true) },
144 ],
145 )}
146 />
147 </>
148 )
149}
151function ClickPanels({ d }: { d: ClickOut }) {
152 return (
153 <>
154 <Row>
155 <SeriesPlot
156 title="the first arrival, after 8.6 m"
157 height={190}
158 xLabel="ms after the exact arrival"
159 yDomain={[-0.6, 1.1]}
160 series={[
161 { x: d.tw1, y: d.ye1, color: ink.reference, label: 'exact', width: 2.4 },
162 { x: d.tw1, y: d.y1, color: method.fd2, label: 'leapfrog', width: 1.3 },
163 ]}
164 />
165 <SeriesPlot
166 title="the twentieth arrival, after 334 m"
167 height={190}
168 xLabel="ms after the exact arrival"
169 yDomain={[-0.6, 1.1]}
170 series={[
171 { x: d.tw20, y: d.ye20, color: ink.reference, label: 'exact', width: 2.4 },
172 { x: d.tw20, y: d.y20, color: method.fd2, label: 'leapfrog', width: 1.3 },
173 ]}
174 />
175 </Row>
176 <Field
177 z={d.spec}
178 xDomain={d.specT}
179 yDomain={d.specF}
180 height={260}
181 kind="seq"
182 reverse
183 title="spectrogram of the leapfrog recording: time in seconds across, frequency in hertz up, level in dB"
184 />
185 <Listen
186 items={[
187 { label: 'exact', clip: d.clipExact },
188 { label: 'leapfrog', clip: d.clipFd },
189 ]}
190 note="the same scale; the exact clip is a plain click train"
191 />
192 </>
193 )
194}
196export const ch02: Chapter = {
197 id: 'leapfrog',
198 nav: 'Leapfrog',
199 part: 'Strings',
200 title: 'Leapfrog, and the magic time step',
201 blurb:
202 'The standard explicit scheme, its stability limit, the one Courant number at which it is exact, and what it sounds like everywhere else.',
203 text,
204 labs: [
205 {
206 id: 'run',
207 title: 'A plucked string at the audio rate',
208 caption:
209 'The time step is one sample at 44.1 kHz, so the grid sets the Courant number: λ = N/200. Try N = 200, then 100 and 50, then 201.',
210 driver,
211 script,
212 controls: [
213 {
214 kind: 'slider',
215 key: 'N',
216 label: 'grid intervals N',
217 min: 0,
218 max: 0,
219 values: [25, 50, 100, 150, 200, 201],
220 format: (v: number) => `${v} (λ = ${(v / 200).toFixed(3)})`,
221 },
222 {
223 kind: 'select',
224 key: 'where',
225 label: 'listen to',
226 options: [
227 { value: 'bridge', label: 'the first node off the bridge' },
228 { value: 'pickup', label: 'a pickup at x = 0.12' },
229 ],
230 },
231 ],
232 defaults: { N: 100, where: 'bridge' },
233 derive: (p: Params) => ({
234 mode: 'run',
235 fs: FS,
236 f0: F0,
237 N: p.N,
238 T: 2.5,
239 where: p.where,
240 beta: 0.2,
241 xp: 0.12,
242 }),
243 editorOpen: true,
244 render: (d: RunOut) => <RunPanels d={d} />,
245 summary: (d: RunOut) =>
246 d.unstable
247 ? [
248 { label: 'Courant number λ', value: fmtNum(d.lam, 4) },
249 { label: 'blew up after', value: `${fmtNum(d.blowup, 3)} ms` },
250 ]
251 : [
252 { label: 'Courant number λ', value: fmtNum(d.lam, 4) },
253 { label: 'rms error, last second', value: fmtErr(d.relerr) },
254 { label: 'partial 1, predicted', value: fmtCents(d.pred1) },
255 { label: 'partial 10, predicted', value: fmtCents(d.pred10) },
256 ],
257 },
258 {
259 id: 'cents',
260 title: 'The pitch of every partial, measured',
261 caption:
262 'Each marker is a frequency read off the simulation; each line is the dispersion relation. Change the grid and watch the lower plot not move.',
263 driver,
264 script,
265 controls: [{ kind: 'slider', key: 'N', label: 'grid intervals N', min: 0, max: 0, values: [25, 50, 100, 200] }],
266 defaults: { N: 100 },
267 derive: (p: Params) => ({ mode: 'cents', f0: F0, N: p.N, lams: [0.25, 0.5, 0.75, 0.9, 1] }),
268 render: (d: CentsOut) => <CentsPanels d={d} />,
269 summary: (d: CentsOut) => [
270 { label: 'measured minus predicted, largest', value: fmtCents(d.agree) },
271 { label: 'error at λ = 1, largest', value: fmtCents(Math.max(...d.meas[d.meas.length - 1].map((v) => Math.abs(v ?? 0)))) },
272 ],
273 },
274 {
275 id: 'click',
276 title: 'The sound of numerical dispersion',
277 caption:
278 'A click travels round a 17.15 m loop of tube and passes a microphone twenty times a second. At λ = 1 it stays a click; below it, listen.',
279 driver,
280 script,
281 controls: [
282 {
283 kind: 'select',
284 key: 'lam',
285 label: 'Courant number λ',
286 options: [
287 { value: 1, label: '1 (the magic time step)' },
288 { value: 0.75, label: '0.75' },
289 { value: 0.5, label: '0.5' },
290 { value: 0.25, label: '0.25' },
291 ],
292 },
293 ],
294 defaults: { lam: 0.5 },
295 derive: (p: Params) => ({ mode: 'click', fs: FS, c: 343, L: 17.15, lam: p.lam, T: 1.2, sigt: 3e-5 }),
296 render: (d: ClickOut) => <ClickPanels d={d} />,
297 summary: (d: ClickOut) => [
298 { label: 'grid intervals', value: String(d.N) },
299 { label: 'highest frequency the grid carries', value: `${(d.fcut / 1000).toFixed(2)} kHz` },
300 { label: 'relative rms error', value: fmtErr(d.rel) },
301 ],
302 },
303 ],
304}