Plate: Lyapunov exponent (shadow plate; real-time chart in a new Real time card, note chart and λ vs E off by default); Parameters folded under the plate
7 changed files+207−34
src/plate/README.mdmodified+4−2View file
@@ -1,8 +1,8 @@
11 # plate
22
3-Struck square plate, step 4 of "Listening for chaos" (after double-pendulum, fput-chain and beam), the first two-dimensional system on the way to cymbal and gong. Work in progress: step 1 (2026-10-09), a linear, simply supported square plate with mallet strike, animation and sped-up note with spectrogram; step 2 (2026-10-09), von Kármán stretching; step 3 (2026-10-09), free edges; step 4 (2026-10-09), stretching with free edges (and the free layout moved to the grid vertices); step 5 (2026-10-09), clamped edges.
3+Struck square plate, step 4 of "Listening for chaos" (after double-pendulum, fput-chain and beam), the first two-dimensional system on the way to cymbal and gong. Work in progress: step 1 (2026-10-09), a linear, simply supported square plate with mallet strike, animation and sped-up note with spectrogram; step 2 (2026-10-09), von Kármán stretching; step 3 (2026-10-09), free edges; step 4 (2026-10-09), stretching with free edges (and the free layout moved to the grid vertices); step 5 (2026-10-09), clamped edges; step 6 (2026-10-09), Lyapunov exponent.
44
5-Current state: N × N masses (N per side, 2 to 32, default 8) on vertical rails at the nodes of a square grid, spacing a, moving up and down only. Bending κ on the five-point Laplacian at each node, energy κ (L u)² / 2, force −κ L(L u). The edge nodes are supports at u = 0 with no bending term (mirror ghost), which makes the edge simply supported. Stretching (von Kármán), below, set by the thickness h (default 5 cm; 0 = linear). Damping as in the beam: drag (τ, all modes alike) and a frequency-dependent loss (dashpots alongside the grid springs, force −c L v, set by τ top at the top of the band). Parameters are given for the 8 × 8 plate (m 1 kg, κ 20 N/m, h 0.05 m, a 0.5 m, so the side is 4.5 m; τ 540 s, τ top 200 s); Masses per side refines the same plate (`refine`: m with the cell area, κ inversely, the loss coefficient c unchanged). Strike: a kick (pyramid velocity bump of half-width 0.5 m) or the shared Hertz mallet (`../common/mallet.ts`, M 1 kg, K 1000) over the same bump, at (x, y) m from the lower left corner (default (1.3, 1.1)); clicking the plate moves the strike point. The card shows the plate from above as a height map (red up, blue down, scaled to the recent largest |u|), with the strike point (cross) and the listening point (circle). The sound panel is the beam's, reduced: velocity, displacement or acceleration at a listening point (default (3.1, 1.7), bilinear between masses), pitch of the lowest mode, the whole decay as a note through the shared renderer and anti-alias filter, spectrogram with dashed lines at the distinct mode pitches, energy chart, Replay and Save WAV. Not yet: Lyapunov exponent, ghost, mode bars and histories, cross-section.
5+Current state: N × N masses (N per side, 2 to 32, default 8) on vertical rails at the nodes of a square grid, spacing a, moving up and down only. Bending κ on the five-point Laplacian at each node, energy κ (L u)² / 2, force −κ L(L u). The edge nodes are supports at u = 0 with no bending term (mirror ghost), which makes the edge simply supported. Stretching (von Kármán), below, set by the thickness h (default 5 cm; 0 = linear). Damping as in the beam: drag (τ, all modes alike) and a frequency-dependent loss (dashpots alongside the grid springs, force −c L v, set by τ top at the top of the band). Parameters are given for the 8 × 8 plate (m 1 kg, κ 20 N/m, h 0.05 m, a 0.5 m, so the side is 4.5 m; τ 540 s, τ top 200 s); Masses per side refines the same plate (`refine`: m with the cell area, κ inversely, the loss coefficient c unchanged). Strike: a kick (pyramid velocity bump of half-width 0.5 m) or the shared Hertz mallet (`../common/mallet.ts`, M 1 kg, K 1000) over the same bump, at (x, y) m from the lower left corner (default (1.3, 1.1)); clicking the plate moves the strike point. The card shows the plate from above as a height map (red up, blue down, scaled to the recent largest |u|), with the strike point (cross) and the listening point (circle). The sound panel is the beam's, reduced: velocity, displacement or acceleration at a listening point (default (3.1, 1.7), bilinear between masses), pitch of the lowest mode, the whole decay as a note through the shared renderer and anti-alias filter, spectrogram with dashed lines at the distinct mode pitches, energy chart, Replay and Save WAV. Not yet: ghost, mode bars and histories, cross-section.
66
77 Modes: L is diagonal in the discrete sine basis, so the modes are exact, S_mn ∝ sin(mπ(i + 1)/(N + 1)) sin(nπ(j + 1)/(N + 1)), with ω_mn = √(κ/m) (l_m + l_n), l_m = 4 sin²(mπ/(2(N + 1))) (`Modes`, sorted by frequency with (m, n)). For low modes ω_mn ∝ m² + n², the continuum plate; (m, n) and (n, m) are degenerate. Ratios to the lowest, against the continuum's 2.5 and 4.0: (1, 2) 2.309 and (2, 2) 3.618 at N = 4, 2.440 and 3.879 at N = 8, 2.483 and 3.966 at N = 16, 2.496 and 3.991 at N = 32; ω_11 stays at 1.054 to 1.089 rad/s across N. Check (2026-10-09, N = 8, kick at (1.3, 1.1), undamped, 4000 s, displacement at (3.1, 1.7)): 32 of the 33 distinct frequencies appear as spectral peaks within 0.5 cents of the Verlet-corrected prediction (2/h) asin(ωh/2) (FFT resolution 2.5 cents); the missing one, (7, 6), is barely excited there (peak 1e-4 of the largest). Energy wobble 1.2e-3 (bounded). Measured decay times match `decayTimes` (exact, L is diagonal in the modes) to the energy wobble: (1, 1) 498 s, (3, 5) 249 s, (8, 8) 146 s.
88
@@ -25,3 +25,5 @@ Clamped edges (step 5, 2026-10-09). Edges select: simply supported, clamped, fre
2525 A membrane tension T for held edges (drum-like: springs between neighbours and to the supports, energy T/2 uᵀLu, exact sine modes with ω² = (T l + κ l²)/m when simply supported) was added in step 5 and removed the same day: the suggested drum setting (κ 2, T 200, h 1 mm) did not sound like a drum. Drums to be revisited later.
2626
2727 Known cost issue (2026-10-09): the 4x oversampled note forces at least 4 steps per audio sample, and the step stability needs is longer than that whenever the speed-up is small, that is when the lowest mode is close in frequency to the top of the band: clamped (lowest mode higher; a tensioned membrane, tried briefly, was far worse, about 8x more steps than needed). Clamped default note: 12 s to render in the browser. Fix planned: read the signal between steps by interpolation (shared note renderer).
28+
29+Lyapunov exponent (step 6, 2026-10-09): as the beam's (two-trajectory method, `../common/lyapunov.ts`), with the plate's shadow in `lyapunov.ts`: a second plate D0 = 1e-8 away, distances with displacements in units of A/√(N²) (A the lowest mode's amplitude at which stretching equals bending, `nonlinearAmplitude`; 1 m for a linear plate) and momenta in m ω₁ times that; the initial offset is a fixed pseudo-random pattern (free: its rigid part removed). Layout now as the beam's: Parameters folded under the plate, a Real time card on the right with Lyapunov (λ 10 s and since the strike, first 5 s left out of the y range) and Energy error charts. Sound panel: Lyapunov chart and λ vs E scatter, both off at first (the shadow doubles the time to compute a note: default note 4.3 s without, about 12 s with, in the browser). Exponents (2026-10-09, N = 8, kick at (1.3, 1.1), undamped, 300 s, λ since the strike, 1/s): linear plate 0.007 at every energy (finite-time; the last 100 s, −0.001: zero); simply supported, h 5 cm: 0.008, 0.065, 0.26, 0.53, 1.01 at 0.001, 0.1, 1, 10, 100 J; free, h 5 cm: 0.006, 0.015, 0.16, 1.54, 8.4. The free plate's exponent grows much faster with energy.
src/plate/index.htmlmodified+31−17View file
@@ -50,24 +50,29 @@
5050 </select>
5151 </div>
5252 </div>
53+ <details>
54+ <summary>Parameters <button type="button" class="help" data-help="params" title="About the parameters" aria-label="About the parameters">?</button></summary>
55+ <div id="param-form">
56+ <label>m <input id="p-m" type="number" step="any" /> <span class="unit">kg</span></label>
57+ <label>κ <input id="p-kappa" type="number" step="any" /> <span class="unit">N/m (bending)</span></label>
58+ <label>h <input id="p-thick" type="number" step="any" /> <span class="unit">m (thickness, sets the stretching; 0 = linear)</span></label>
59+ <label>a <input id="p-a" type="number" step="any" /> <span class="unit">m (spacing)</span></label>
60+ <label>τ <input id="p-tau" type="number" step="any" /> <span class="unit">s (decay, 0 = off)</span></label>
61+ <label>τ top <input id="p-tautop" type="number" step="any" /> <span class="unit">s (decay of the highest modes from the loss, 0 = off)</span></label>
62+ <label>mallet M <input id="p-mallet-m" type="number" step="any" /> <span class="unit">kg</span></label>
63+ <label>mallet K <input id="p-mallet-k" type="number" step="any" /> <span class="unit">N/m^1.5 (contact stiffness)</span></label>
64+ </div>
65+ <div id="param-actions" class="row">
66+ <button id="param-apply">Apply</button>
67+ <span id="param-status"></span>
68+ </div>
69+ <div id="params"></div>
70+ </details>
5371 </section>
5472 <section class="card">
55- <div class="card-head">Parameters <button type="button" class="help" data-help="params" title="About the parameters" aria-label="About the parameters">?</button></div>
56- <div id="param-form">
57- <label>m <input id="p-m" type="number" step="any" /> <span class="unit">kg</span></label>
58- <label>κ <input id="p-kappa" type="number" step="any" /> <span class="unit">N/m (bending)</span></label>
59- <label>h <input id="p-thick" type="number" step="any" /> <span class="unit">m (thickness, sets the stretching; 0 = linear)</span></label>
60- <label>a <input id="p-a" type="number" step="any" /> <span class="unit">m (spacing)</span></label>
61- <label>τ <input id="p-tau" type="number" step="any" /> <span class="unit">s (decay, 0 = off)</span></label>
62- <label>τ top <input id="p-tautop" type="number" step="any" /> <span class="unit">s (decay of the highest modes from the loss, 0 = off)</span></label>
63- <label>mallet M <input id="p-mallet-m" type="number" step="any" /> <span class="unit">kg</span></label>
64- <label>mallet K <input id="p-mallet-k" type="number" step="any" /> <span class="unit">N/m^1.5 (contact stiffness)</span></label>
65- </div>
66- <div id="param-actions" class="row">
67- <button id="param-apply">Apply</button>
68- <span id="param-status"></span>
69- </div>
70- <div id="params"></div>
73+ <div class="card-head">Real time <button type="button" class="help" data-help="realtime" title="About the real-time charts" aria-label="About the real-time charts">?</button></div>
74+ <div class="chips"><span id="toggles"></span></div>
75+ <div id="traces"></div>
7176 </section>
7277 </div>
7378 <details class="card panel" id="sound" open>
@@ -94,6 +99,10 @@
9499 <p id="sound-hint" class="muted">Strike to hear the whole decay.</p>
95100 <figure id="spec-fig"><div id="spectrogram"></div></figure>
96101 <div id="sound-traces"></div>
102+ <figure id="scatter-fig" hidden>
103+ <div id="scatter"></div>
104+ <div class="row scatter-actions"><button id="scatter-clear" type="button">Clear</button></div>
105+ </figure>
97106 </details>
98107 </main>
99108 <dialog id="help">
@@ -113,10 +122,15 @@
113122 <p>The stretching follows the von Kármán plate equations: the in-plane stress comes from a stress function, solved on the grid at every step. A simply supported edge may move outward but not along itself; a free edge is free in its own plane too. The nonlinear scale below is the energy of the lowest mode at which its stretching energy equals its bending energy.</p>
114123 <p>With free edges the bending energy also includes the plate's twisting term (Poisson's ratio 0.3), which sets the edge conditions; the modes are computed, and their frequencies converge to plate theory (Leissa's tables) as the spacing squared (the lowest is 2.6% low at 9 × 9, 0.2% at 33 × 33). On a free plate the outer masses sit on the edge, carrying half a cell's mass (a quarter at the corners), so the free plate's side is one spacing less than the number of masses. Clamped edges hold the slope at zero too; the modes are computed (no (m, n)), and the plate is stiffer than a simply supported one, with overtones about as widely spread (plate theory: 1, 2.04, 3.01, 3.66, ...). For a simply supported plate, the normal modes are products of sines, numbered (m, n) by their half waves along x and y. Their frequencies go as m² + n² for the low ones (the continuum plate), so (1, 2) and (2, 1) share a pitch, and the ratios to the lowest, 2.5, 4, 5, 6.5, 8.5, ..., are not harmonic. The derived values list the lowest modes and how close the grid is to the continuum.</p>
115124 </template>
125+ <template id="help-realtime">
126+ <h2>Real time</h2>
127+ <p>Charts of the plate on the left, on the same clock. Click a chip to show or hide one; hover to read values, click to pin the cursor.</p>
128+ <p>Lyapunov is the rate (1/s) at which nearby motions separate, measured with a hidden second plate kept a tiny distance away: positive for chaos, near 0 for regular motion (a linear plate, or a weak strike). Shown over the last 10 s and since the strike; the first 5 s are left out of the y range. Energy error should stay tiny.</p>
129+ </template>
116130 <template id="help-sound">
117131 <h2>Sped up · listen</h2>
118132 <p>A copy of the same strike, run fast enough to be heard: Pitch sets where the lowest normal mode lands, and the whole decay plays in a few seconds. While this panel is open each Strike plays a note; Replay plays the last one again with the current settings.</p>
119- <p>Signal chooses what is heard, at a listening point (x, y in metres, between masses by interpolation). Save WAV saves the last note as heard (16-bit mono, normalized). The spectrogram has dashed lines at the normal-mode pitches. Notes are anti-alias filtered (motion above 20 kHz at 48 kHz is removed rather than folded back); a note says when the top modes are above that and so are not heard.</p>
133+ <p>Signal chooses what is heard, at a listening point (x, y in metres, between masses by interpolation). Save WAV saves the last note as heard (16-bit mono, normalized). The spectrogram has dashed lines at the normal-mode pitches. The Lyapunov chip shows the exponent over the last 10 s of the plate's own time through the note, and λ vs E collects it against the energy over notes, earlier notes faded (Clear starts over); both need a hidden second plate run along, which doubles the time to compute a note, so they are off at first. Notes are anti-alias filtered (motion above 20 kHz at 48 kHz is removed rather than folded back); a note says when the top modes are above that and so are not heard.</p>
120134 </template>
121135 <script type="module" src="./main.ts"></script>
122136 </body>
src/plate/lyapunov.tsadded+63−0View file
@@ -0,0 +1,63 @@
1+// The plate's Lyapunov shadow (see ../common/lyapunov.ts): a second plate
2+// kept D0 from the plate.
3+//
4+// Distance: displacements in units of ur = A / sqrt(N^2), with A the
5+// lowest mode's amplitude at which stretching matches its bending energy
6+// (the size of a nonlinear motion, per mass), and momenta in units of m
7+// omega_1 ur (the lowest mode's momentum at that size), so both count
8+// comparably. A linear plate (no stretching) has no such size; ur = 1 m
9+// then (its exponent is 0 whatever the units).
10+
11+import { D0, type Shadow } from "../common/lyapunov";
12+import { Verlet, copyState, nonlinearAmplitude, removeRigid, type Modes, type Params, type State } from "./plate";
13+
14+// Fixed pseudo-random pattern in [-1, 1] (so every run starts the shadow
15+// along the same direction).
16+function pattern(n: number, seed: number): Float64Array {
17+ const v = new Float64Array(n);
18+ let x = seed >>> 0;
19+ for (let k = 0; k < n; k++) {
20+ x = (Math.imul(x, 1664525) + 1013904223) >>> 0;
21+ v[k] = (x / 2 ** 32) * 2 - 1;
22+ }
23+ return v;
24+}
25+
26+// A shadow of the plate's state s, for Lyapunov.start. s is read at every
27+// step, so it must stay the plate's state object until the next start.
28+export function shadowOf(P: Params, modes: Modes, s: State): Shadow {
29+ const n = s.u.length, A = nonlinearAmplitude(P, modes);
30+ const ur = Number.isFinite(A) ? A / Math.sqrt(n) : 1, pr = P.m * modes.omega[0] * ur;
31+ // Initial separation: a fixed pattern in u and in p (free: without its
32+ // rigid part, which the plate does not have), scaled to D0.
33+ const du = pattern(n, 1), dp = pattern(n, 2);
34+ if (P.edges === "free") {
35+ removeRigid(P, du);
36+ removeRigid(P, dp, true);
37+ }
38+ let d2 = 0;
39+ for (let k = 0; k < n; k++) d2 += du[k] * du[k] + dp[k] * dp[k];
40+ const f = D0 / Math.sqrt(d2);
41+ const sh = copyState(s);
42+ for (let k = 0; k < n; k++) {
43+ sh.u[k] += f * du[k] * ur;
44+ sh.p[k] += f * dp[k] * pr;
45+ }
46+ const verlet = new Verlet(P, sh);
47+
48+ return {
49+ step: (h) => verlet.step(sh, h),
50+ separation: () => {
51+ let d = 0;
52+ for (let k = 0; k < n; k++) d += ((sh.u[k] - s.u[k]) / ur) ** 2 + ((sh.p[k] - s.p[k]) / pr) ** 2;
53+ return Math.sqrt(d);
54+ },
55+ pullBack: (r) => {
56+ for (let k = 0; k < n; k++) {
57+ sh.u[k] = s.u[k] + (sh.u[k] - s.u[k]) * r;
58+ sh.p[k] = s.p[k] + (sh.p[k] - s.p[k]) * r;
59+ }
60+ verlet.reset(sh); // forces at the moved displacements
61+ },
62+ };
63+}
src/plate/main.tsmodified+78−5View file
@@ -11,11 +11,14 @@ import { setupHelp } from "../common/help";
1111 import { setupMute } from "../common/mute";
1212 import { PASS } from "../common/antialias";
1313 import { malletStrike, type MalletParams } from "../common/mallet";
14-import { DIAG_EVERY, END_RATIO, UNDAMPED_LENGTH } from "../common/note-render";
14+import { Lyapunov } from "../common/lyapunov";
15+import { DIAG_EVERY, END_RATIO, LYAP_WINDOW, UNDAMPED_LENGTH } from "../common/note-render";
16+import { ScatterPlot } from "../common/scatter";
1517 import { Player } from "../common/sound";
1618 import { NoteSpectrogram } from "../common/spectrogram";
1719 import { TracePanel } from "../common/traces";
1820 import { saveFile, wavBlob } from "../common/wav";
21+import { shadowOf } from "./lyapunov";
1922 import type { RenderRequest, RenderResult } from "./render-worker";
2023 import { SIGNALS, signalByKey } from "./signals";
2124
@@ -128,6 +131,9 @@ document.getElementById("strike")!.addEventListener("click", () => {
128131 showDerived();
129132 t = 0;
130133 ampl = 0;
134+ lyap.start(shadowOf(P, modes, s), t);
135+ traces.clear();
136+ traces.sample(t);
131137 replayButton.disabled = false;
132138 renderNote();
133139 pauseButton.disabled = false;
@@ -148,6 +154,32 @@ plateCanvas.addEventListener("click", (e) => {
148154 // Energy balance: energy plus what damping has dissipated, against the strike's.
149155 const relErr = () => (E0 > 0 ? (energy(P, s) + s.W - E0) / E0 : 0);
150156
157+// Largest Lyapunov exponent of the animated plate, measured from each strike
158+// (two-trajectory method, shadow in lyapunov.ts): over the last 10 s and
159+// since the strike.
160+const lyap = new Lyapunov(LYAP_WINDOW);
161+// Time after a strike left out of the Lyapunov charts' y range: the early
162+// values swing widely before they settle.
163+const LYAP_SKIP = 5;
164+const SAMPLE_DT = 0.01; // s between plotted samples
165+
166+// Real-time charts, sampled every SAMPLE_DT.
167+const traces = new TracePanel(
168+ document.getElementById("toggles")!,
169+ document.getElementById("traces")!,
170+ [
171+ {
172+ label: "Lyapunov", unit: "1/s", window: 60, skip: () => LYAP_SKIP, shown: true,
173+ series: [
174+ { label: "λ 10 s", value: () => lyap.windowed() },
175+ { label: "λ since strike", value: () => lyap.mean() },
176+ ],
177+ },
178+ { label: "Energy error", series: [{ label: "energy balance error", value: relErr }] },
179+ ],
180+ 10,
181+);
182+
151183 // Parameters: numbers in the form; Apply restarts at rest. N is set above.
152184 const paramInputs = {
153185 m: document.getElementById("p-m") as HTMLInputElement,
@@ -263,15 +295,40 @@ const specToggle = document.getElementById("spec-toggle") as HTMLInputElement;
263295 specToggle.addEventListener("change", () => (specFig.hidden = !specToggle.checked));
264296 const specLog = document.getElementById("spec-log") as HTMLInputElement;
265297 specLog.addEventListener("change", () => spec?.setLog(specLog.checked));
266-let noteEnergy = NaN;
298+let noteEnergy = NaN, noteLambda = NaN;
299+// The note's Lyapunov exponent needs a shadow run along in the worker, which
300+// doubles the time to compute a note: off at first.
267301 const soundTraces = new TracePanel(
268302 document.getElementById("sound-chips")!,
269303 document.getElementById("sound-traces")!,
270- [{ label: "Energy", unit: "J", log: true, series: [{ label: "E", value: () => noteEnergy }] }],
304+ [
305+ {
306+ label: "Lyapunov", unit: "1/s",
307+ skip: () => LYAP_SKIP / speedup(), // the same physical time as in the animation, in s of sound
308+ series: [{ label: "λ 10 s", value: () => noteLambda }],
309+ },
310+ { label: "Energy", unit: "J", log: true, series: [{ label: "E", value: () => noteEnergy }] },
311+ ],
271312 10, 960, soundCursor,
272313 );
314+// Lambda against energy, accumulated over notes (chip after the trace chips).
315+const scatter = new ScatterPlot(960);
316+document.getElementById("scatter")!.appendChild(scatter.canvas);
317+const scatterFig = document.getElementById("scatter-fig")!;
318+{
319+ const box = document.createElement("input");
320+ box.type = "checkbox";
321+ box.addEventListener("change", () => (scatterFig.hidden = !box.checked));
322+ const text = document.createElement("span");
323+ text.textContent = "λ vs E";
324+ const chip = document.createElement("label");
325+ chip.className = "chip";
326+ chip.append(box, text);
327+ document.getElementById("sound-chips")!.appendChild(chip);
328+}
329+document.getElementById("scatter-clear")!.addEventListener("click", () => scatter.clear());
273330 let spec: NoteSpectrogram | null = null;
274-let note: { x: Float32Array; energy: Float32Array; rate: number } | null = null;
331+let note: { x: Float32Array; energy: Float32Array; lambda: Float32Array; rate: number } | null = null;
275332 let fed = 0; // samples of the note given to the displays
276333 function feedDisplays(upTo: number) {
277334 if (!note || !spec) return;
@@ -280,7 +337,9 @@ function feedDisplays(upTo: number) {
280337 spec.push(note.x[fed]);
281338 if (fed % DIAG_EVERY === 0) {
282339 noteEnergy = note.energy[fed / DIAG_EVERY];
340+ noteLambda = note.lambda[fed / DIAG_EVERY];
283341 soundTraces.sample(fed / note.rate);
342+ scatter.add(noteEnergy, noteLambda, fed / note.rate);
284343 }
285344 }
286345 }
@@ -305,6 +364,8 @@ function renderNote() {
305364 const req: RenderRequest = {
306365 id: ++requestId, P, state: lastStrike, speedup: speedup(), sampleRate: ctx.sampleRate,
307366 maxLength: 2 * predictedNoteLength(), signal: sig.key, x, y,
367+ // Only when a display needs it (Lyapunov chart or λ vs E).
368+ lyapunov: soundTraces.shownLabels().includes("Lyapunov") || !scatterFig.hidden,
308369 };
309370 noteInfo = `${sig.label} (${x}, ${y}) m`;
310371 noteFile = `plate-${P.edges}-N${P.N}-${strikeLabel}-${energy(P, lastStrike).toPrecision(3)}J-${sig.key}x${x}y${y}m-${pitch().toFixed(0)}Hz.wav`;
@@ -319,7 +380,7 @@ worker.onmessage = (e: MessageEvent<RenderResult>) => {
319380 const rate = player.sampleRate;
320381 const x = e.data.samples;
321382 Player.normalize(x, rate);
322- note = { x, energy: e.data.energy, rate };
383+ note = { x, energy: e.data.energy, lambda: e.data.lambda, rate };
323384 fed = 0;
324385 const span = x.length / rate;
325386 const n = Math.min(16384, 2 ** Math.ceil(Math.log2((2048 * 30) / Math.min(30, pitch()))));
@@ -331,6 +392,7 @@ worker.onmessage = (e: MessageEvent<RenderResult>) => {
331392 document.getElementById("spectrogram")!.replaceChildren(spec.canvas);
332393 soundTraces.clear();
333394 soundTraces.setFixed({ t0: 0, span });
395+ scatter.newNote(span);
334396 const secs = (performance.now() - renderStart) / 1000;
335397 soundStatus.textContent = `${noteInfo} · ${fmtDuration(span)} · computed in ${secs.toFixed(1)} s`;
336398 player.play(x);
@@ -351,6 +413,8 @@ function setParams(next: Params, N: number) {
351413 E0 = 0;
352414 t = 0;
353415 ampl = 0;
416+ lyap.stop();
417+ traces.clear();
354418 nInput.value = String(P.N);
355419 pauseButton.disabled = true;
356420 setRunning(false);
@@ -479,6 +543,7 @@ function drawPlate() {
479543 // Animation: the simulation advances `speed` seconds per second of real time.
480544 let last = performance.now();
481545 let carry = 0;
546+let sinceSample = 0;
482547 let playing = false;
483548 function frame(now: number) {
484549 const dt = Math.min(0.1, (now - last) / 1000);
@@ -487,14 +552,21 @@ function frame(now: number) {
487552 carry += Number(speedSelect.value) * dt;
488553 while (carry >= H) {
489554 verlet.step(s, H);
555+ lyap.step(H, t + H);
490556 carry -= H;
491557 t += H;
558+ sinceSample += H;
559+ if (sinceSample >= SAMPLE_DT) {
560+ sinceSample -= SAMPLE_DT;
561+ traces.sample(t);
562+ }
492563 }
493564 }
494565 let umax = 0;
495566 for (const u of s.u) umax = Math.max(umax, Math.abs(u));
496567 ampl = Math.max(umax, 0.998 * ampl);
497568 drawPlate();
569+ traces.draw();
498570 if (note && soundPanel.open) {
499571 const at = player.position();
500572 if (at !== null) feedDisplays(Math.round(at * note.rate));
@@ -503,6 +575,7 @@ function frame(now: number) {
503575 playing = at !== null;
504576 if (!specFig.hidden) spec?.draw();
505577 soundTraces.draw();
578+ if (!scatterFig.hidden) scatter.draw();
506579 }
507580 requestAnimationFrame(frame);
508581 }
src/plate/plate.tsmodified+21−3View file
@@ -576,6 +576,10 @@ export class Verlet {
576576 this.c = lossCoefficient(P);
577577 forces(P, s.u, this.F, this.tmp);
578578 }
579+ // Call after changing u from outside.
580+ reset(s: State): void {
581+ forces(this.P, s.u, this.F, this.tmp);
582+ }
579583 // Force on each mass at the current displacements.
580584 get force(): Float64Array {
581585 return this.F;
@@ -740,13 +744,27 @@ export function decayTimes(P: Params, modes: Modes): Float64Array {
740744 // (Infinity when linear). Bending goes as the amplitude squared, stretching
741745 // as its fourth power.
742746 export function nonlinearScale(P: Params, modes: Modes): number {
743- if (stretchCoefficient(P) === 0) return Infinity;
747+ const { V2, V4 } = lowestModeEnergies(P, modes);
748+ return V4 > 0 ? (V2 * V2) / V4 : Infinity;
749+}
750+// Bending (V2) and stretching (V4) energy of the lowest mode's shape scaled
751+// to unit length (sum of u^2 = 1 m^2); 0 stretching when linear.
752+function lowestModeEnergies(P: Params, modes: Modes): { V2: number; V4: number } {
744753 const N = P.N;
745754 let u = new Float64Array(N * N);
746755 if (modes.S) u = Float64Array.from(modes.S[0]);
747756 else for (let j = 0; j < N; j++) for (let i = 0; i < N; i++) u[i + N * j] = Math.sin((Math.PI * (i + 1)) / (N + 1)) * Math.sin((Math.PI * (j + 1)) / (N + 1));
748- const V2 = linearEnergy(P, u), V4 = stretching(P, u);
749- return V4 > 0 ? (V2 * V2) / V4 : Infinity;
757+ let n2 = 0;
758+ for (const x of u) n2 += x * x;
759+ for (let k = 0; k < u.length; k++) u[k] /= Math.sqrt(n2);
760+ return { V2: linearEnergy(P, u), V4: stretching(P, u) };
761+}
762+// Nonlinear amplitude (m): the length of the lowest mode's displacement
763+// vector at which its stretching energy equals its bending energy (Infinity
764+// when linear).
765+export function nonlinearAmplitude(P: Params, modes: Modes): number {
766+ const { V2, V4 } = lowestModeEnergies(P, modes);
767+ return V4 > 0 ? Math.sqrt(V2 / V4) : Infinity;
750768 }
751769
752770 // Integrator step (s) for motion with energy E: STEP / omega, omega at the
src/plate/render-worker.tsmodified+6−6View file
@@ -1,9 +1,10 @@
11 // Renders a note: the plate from a given state, sped up, sampled as audio
2-// (the loop and energy: ../common/note-render.ts), the step re-chosen as the
3-// energy falls. No Lyapunov shadow yet.
2+// (the loop, energy and Lyapunov exponent: ../common/note-render.ts), the
3+// step re-chosen as the energy falls.
44
55 import { Modes, Verlet, copyState, energy, isDamped, stepFor, type Params, type State } from "./plate";
66 import { renderNote, type Note } from "../common/note-render";
7+import { shadowOf } from "./lyapunov";
78 import { signalAt, signalByKey } from "./signals";
89
910 export interface RenderRequest {
@@ -16,6 +17,7 @@ export interface RenderRequest {
1617 signal: string; // signal key
1718 x: number; // m, listening point
1819 y: number;
20+ lyapunov: boolean; // run the Lyapunov shadow
1921 }
2022
2123 export interface RenderResult extends Note {
@@ -35,11 +37,9 @@ self.onmessage = (e: MessageEvent<RenderRequest>) => {
3537 stepFor: (E) => stepFor(P, modes, E),
3638 step: (h) => verlet.step(s, h),
3739 signal: () => signalAt(signal, P, s, verlet.force, r.x, r.y),
38- shadow: () => {
39- throw new Error("no Lyapunov shadow for the plate yet");
40- },
40+ shadow: () => shadowOf(P, modes, s),
4141 },
42- { speedup: r.speedup, sampleRate: r.sampleRate, damped: isDamped(P), maxLength: r.maxLength, lyapunov: false },
42+ { speedup: r.speedup, sampleRate: r.sampleRate, damped: isDamped(P), maxLength: r.maxLength, lyapunov: r.lyapunov },
4343 );
4444 const result: RenderResult = { id: r.id, ...note };
4545 (self as unknown as Worker).postMessage(result, [result.samples.buffer, result.energy.buffer, result.lambda.buffer]);
src/plate/style.cssmodified+4−1View file
@@ -6,7 +6,10 @@ canvas { display: block; width: 100%; }
66 [hidden] { display: none !important; }
77 #plate { aspect-ratio: 1; cursor: crosshair; }
88 .card-head { font-size: 12px; margin-bottom: 12px; display: flex; align-items: center; }
9-#sound-traces canvas { height: auto; margin-top: 12px; }
9+#traces canvas, #sound-traces canvas { height: auto; margin-top: 12px; }
10+#toggles { display: contents; }
11+#scatter canvas { width: 100%; height: auto; }
12+#param-form { margin-top: 10px; }
1013 #sound-chips { margin-left: auto; }
1114 input[type="number"] { width: 4.5em; font-family: var(--mono); }
1215 #energy-value { min-width: 4.5em; text-align: right; }