Beam: frequency-dependent loss (dashpots alongside the springs, set by the top-of-band decay time τ top, default off); per-mode decay times in the derived values
5 changed files+101−15
src/beam/README.mdmodified+3−2View file
@@ -1,8 +1,8 @@
11 # beam
22
3-Struck beam with bending stiffness, step 3 of "Listening for chaos" (after double-pendulum and fput-chain). Work in progress: step 1 (physics, animation, sped-up note with spectrogram), step 2 (Lyapunov exponent), step 3 (mode bars and mode histories), step 4 (ghost), step 5 (λ vs E), curvature (rise h of a curved rest shape).
3+Struck beam with bending stiffness, step 3 of "Listening for chaos" (after double-pendulum and fput-chain). Work in progress: step 1 (physics, animation, sped-up note with spectrogram), step 2 (Lyapunov exponent), step 3 (mode bars and mode histories), step 4 (ghost), step 5 (λ vs E), curvature (rise h of a curved rest shape), frequency-dependent loss (τ top).
44
5-Current state: N masses (2 to 64, default 8) on vertical rails, moving sideways only, with tension k, bending κ and a stretching nonlinearity k3 (below). Ends (a select next to Masses; changing it returns the beam to rest): clamped–clamped (default), pinned–pinned, or cantilever (clamped at the left, free at the right). Strike gives energy E (log slider, 10 µJ to 100 J, default 2.51 J) as velocity, on one mass alone (default mass 2) or in the shape of one normal mode. The card draws the beam from the side, at rest shape plus displacement, scaled to fit (scale given; the rest shape is a faint line), with rails as faint lines, a hatched wall for a clamp and a triangle for a pin, plus the time, energy and relative energy error. The Ghost chip (as in fput-chain; controls in `../common/ghost.ts`) adds a second beam started from the same state with every mass moved by 10^k / sqrt(N) m (total offset 10^k m, k = -12 to -1, default -10), run in lockstep without renormalization and drawn in translucent orange underneath; restarted on Strike, on "Restart ghost", when switched on or when the offset changes; a readout gives the separation (root sum of squared displacement differences, m) and the time since the ghost started. Default strike: the separation stays near 1e-10 m for ~10 s, then grows to 1.7e-8 m by 39.5 s, e^5.1, matching λ since the strike (0.13 /s × 39.5 s). A "Real time" card on the right has chips: Mode bars (on; the current harmonic energy of each mode relative to the strike, log scale over 12 decades, with the number of modes sharing the energy, 1 / sum (E_n / sum E)^2, and the mode frequency range; `../common/modebars.ts`), Mode history (off; each mode's share of the harmonic energy over the last 300 s as a heat map, mode 1 at the bottom; `LiveModeHistory` in `../common/modehistory.ts`), and charts (`../common/traces.ts`): Lyapunov (on; the largest Lyapunov exponent over the last 10 s and since the strike, last 60 s shown, the first 5 s after a strike left out of the y range) Mode energies (off; log, last 60 s, values below 1e-12 of the strike left out), Energy error (off) and Separation (log, last 60 s; offered only with the ghost on). The exponent is the two-trajectory method of `../common/lyapunov.ts` with the beam's shadow in `lyapunov.ts`: a hidden second beam 1e-8 away, displacements in units of A / sqrt(N) and momenta in units of m omega_1 A / sqrt(N), A the mode-1 amplitude at the nonlinear scale. Speed 0.1x to 20x real time; Pause or space. Parameters (m, k, κ, k3, rise h, a, τ, stretching form) fold out with derived values (mode frequencies, ratios to mode 1, the nonlinear scale, integrator step); Apply returns the beam to rest.
5+Current state: N masses (2 to 64, default 8) on vertical rails, moving sideways only, with tension k, bending κ and a stretching nonlinearity k3 (below). Ends (a select next to Masses; changing it returns the beam to rest): clamped–clamped (default), pinned–pinned, or cantilever (clamped at the left, free at the right). Strike gives energy E (log slider, 10 µJ to 100 J, default 2.51 J) as velocity, on one mass alone (default mass 2) or in the shape of one normal mode. The card draws the beam from the side, at rest shape plus displacement, scaled to fit (scale given; the rest shape is a faint line), with rails as faint lines, a hatched wall for a clamp and a triangle for a pin, plus the time, energy and relative energy error. The Ghost chip (as in fput-chain; controls in `../common/ghost.ts`) adds a second beam started from the same state with every mass moved by 10^k / sqrt(N) m (total offset 10^k m, k = -12 to -1, default -10), run in lockstep without renormalization and drawn in translucent orange underneath; restarted on Strike, on "Restart ghost", when switched on or when the offset changes; a readout gives the separation (root sum of squared displacement differences, m) and the time since the ghost started. Default strike: the separation stays near 1e-10 m for ~10 s, then grows to 1.7e-8 m by 39.5 s, e^5.1, matching λ since the strike (0.13 /s × 39.5 s). A "Real time" card on the right has chips: Mode bars (on; the current harmonic energy of each mode relative to the strike, log scale over 12 decades, with the number of modes sharing the energy, 1 / sum (E_n / sum E)^2, and the mode frequency range; `../common/modebars.ts`), Mode history (off; each mode's share of the harmonic energy over the last 300 s as a heat map, mode 1 at the bottom; `LiveModeHistory` in `../common/modehistory.ts`), and charts (`../common/traces.ts`): Lyapunov (on; the largest Lyapunov exponent over the last 10 s and since the strike, last 60 s shown, the first 5 s after a strike left out of the y range) Mode energies (off; log, last 60 s, values below 1e-12 of the strike left out), Energy error (off) and Separation (log, last 60 s; offered only with the ghost on). The exponent is the two-trajectory method of `../common/lyapunov.ts` with the beam's shadow in `lyapunov.ts`: a hidden second beam 1e-8 away, displacements in units of A / sqrt(N) and momenta in units of m omega_1 A / sqrt(N), A the mode-1 amplitude at the nonlinear scale. Speed 0.1x to 20x real time; Pause or space. Parameters (m, k, κ, k3, rise h, a, τ, τ top, stretching form) fold out with derived values (mode frequencies, ratios to mode 1, the nonlinear scale, integrator step); Apply returns the beam to rest.
66
77 "Sped up · listen" (as in fput-chain): each Strike renders the same strike in a worker (`render-worker.ts`), sped up so the lowest mode sounds at the Pitch (5 to 1000 Hz, default 100), and plays it; the page starts muted. A damped note is the whole decay to -60 dB. A red warning shows when the top mode is above half the sample rate (it would alias). Signals (`signals.ts`): velocity (default), displacement or acceleration of mass j (default 2), or the sum of accelerations (the net force of the supports). Replay, Save WAV (named e.g. beam-clamped-local-N8-mass2-2.51J-v2-100Hz.wav), spectrogram over the whole note with dashed lines at the mode pitches. Chips: Spectrogram (on; frequency axis logarithmic by default, from a round number at or below half the lowest mode pitch, with a "log f" chip to switch to linear from 0; each pixel row shows the loudest FFT bin in its band; mode lines each drawn only while at least 6 px from their neighbours: at the default note the six lowest are drawn, modes 7 and 8 are 2 px apart; on the old linear axis to 10 kHz the modes filled the bottom eighth), Modes (off; the note's mode-energy shares over the whole note, one column per 256 samples, from the worker; `ModeHistory` in `../common/modehistory.ts`), Lyapunov (on; the note's exponent over the last 10 s of physical time, from a shadow run in the worker; the first 5 s of physical time left out of the y range; the shadow runs only while this chart or λ vs E is shown, as it doubles the cost) Energy (log) and λ vs E (off; `../common/scatter.ts` as in fput-chain: one point per 256 samples, accumulated over notes, log energy axis, colored by time into the note, earlier notes faded, hover for values, Clear; not cleared when parameters change). All but λ vs E span the whole note and share one cursor (playhead while playing; hover, click to pin). The render loop is shared with fput-chain (`../common/note-render.ts`). The status line gives the render time (default note 1.0 to 1.4 s with the shadow, 0.6 s without).
88
@@ -20,6 +20,7 @@ Live site: https://chaos-apps--listening-for-chaos.vault1-sites.magland.org/beam
2020 - Normal modes of the linear beam from the stiffness matrix (assembled from the force routine) by Jacobi; mode energies as in fput-chain. Checked: pinned frequencies equal the closed form omega_n^2 = (4 k s^2 + 16 κ s^4) / m, s = sin(pi n / (2 (N + 1))), to 1e-14; N = 32 ratios clamped 2.75, 5.35, 8.78 (beam theory 2.76, 5.40, 8.93), cantilever 6.26, 17.48, 34.09 (6.27, 17.55, 34.39). At N = 8: clamped 1, 2.59, 4.68, 6.97, top 13.2; pinned 1, 3.88, 8.29, top 32.2; cantilever 1, 6.13, 16.4, top 75.
2121 - Nonlinear scale: the energy of mode 1 alone at the amplitude where its stretching energy equals its linear energy (N = 8, local: clamped 176 mJ, pinned 8.7 mJ, cantilever 1.1 mJ).
2222 - Damping: drag proportional to mass, every mode's amplitude decays as exp(-t / τ); applied exactly in half steps around each Verlet step, dissipated energy in W.
23+- Frequency-dependent loss (default off): a dashpot alongside each spring, force -c (L v) with L the springs' second-difference matrix (the discrete form of Bilbao's plate loss term), so short waves lose more. Set by τ top, the amplitude decay time it gives the top of the band, c = 2 m / (τ top λmax(L)); a mode with q = S L S decays at c q / 2m on top of the drag. The derived values list each mode's decay time (drag and loss together; for clamped and cantilever ends the loss is treated as diagonal in the modes). Integrated as two Crank-Nicolson half steps around the Verlet step (tridiagonal, O(N), never adds energy); the dissipated energy goes into W. Checked (Node): τ top = 0 bit-identical to before (all ends, both forms, straight and curved); predicted and measured mode decay times agree to 1 to 2% (N = 8, τ top 5 s: clamped 122 s for mode 1 down to 5 s for mode 8); energy balance <= 1.5e-5 (N = 16, 10 J, 300 T1); N = 64 step 1.0 -> 2.4 µs. Default strike with τ top = 50 s (decay times 374 s for mode 1 down to 46 s): the upper lines die first and the note settles onto mode 1; spectral centroid 1000 Hz at 0.1 s, 160 Hz at 1 s, 100 Hz from 2 s (without loss 860 Hz at 1 s, 840 Hz at 3 s); note 3.6 s instead of 7.8 s. The header's predicted note length uses the slowest mode, so it overestimates when most of the energy starts in fast-decaying modes (5.5 s predicted here).
2324 - Integrator: velocity Verlet, step 0.03 / omega_eff, omega_eff^2 = (Kmax + 12 k3 d^2) / m, Kmax the top eigenvalue of the stiffness matrix and d = (4 E / k3)^(1/4) (Kirchhoff: the k3 term over sqrt(nS)). Relative energy error <= 8e-5 over 200 mode-1 periods, all ends and both forms, 10 mJ to 100 J.
2425 - Default note (clamped, local, 2.51 J on mass 2, τ 540 s, 100 Hz): λ ~0.27 /s at the strike (λ T1 ~1.3), falling into the noise of the 10 s window (±0.05 /s) by ~1.6 s of sound (E ~0.14 J), as the spectrogram resolves into lines. Kirchhoff at the same strike: ~0.04 /s, then noise. A 50 mJ strike: noise throughout. Cantilever: ~0.2 /s for the first ~0.3 s of its 1.3 s note. On Modes, the shares shuffle while λ is up and freeze into steady stripes as it falls (~2 s of sound); on Mode bars about 6 to 7 of the 8 modes share the energy early on. λ vs E over a 2.51 J and a 25 J note (clamped, local): λ in the noise around 0 at low energy, above it from roughly 0.3 J, ~1.4 /s at 25 J (25 J note: 7.5 s, computed in 2.1 s).
2526 - First look (offline, N = 8, hammer on mass 2, no damping, two-trajectory Lyapunov exponent times the mode-1 period T1): local stretching is chaotic from ~1 J (clamped 0.57 at 1 J, 3.3 at 10 J; pinned 1.8 at 1 J; cantilever 1.2 at 0.1 J). Kirchhoff is much less chaotic at the same energy (clamped 0.04 at 1 J, 0.24 at 10 J; pinned ~0 up to 1 J, 0.09 at 10 J).
src/beam/beam.tsmodified+75−1View file
@@ -29,6 +29,14 @@
2929 // Optional damping: linear drag on each mass proportional to its mass,
3030 // force -gamma m v = -gamma p, so every linear mode's amplitude decays as
3131 // exp(-t / tau) with tau = 2 / gamma.
32+// Optional frequency-dependent loss: a dashpot alongside each spring, force
33+// -c (v[j] - v[j - 1]) on its ends, i.e. -c (L v) with L the second-difference
34+// (graph Laplacian) matrix of the springs, as in Bilbao's plate loss term.
35+// It damps short waves more: a mode with q = S L S (its shape S) loses
36+// amplitude at rate c q / 2m; q grows with the wavenumber squared, so for a
37+// beam (omega ~ wavenumber^2) the rate grows roughly with frequency.
38+// Set by tauTop, the amplitude decay time it gives the top of the band
39+// (q = the largest eigenvalue of L): c = 2 m / (tauTop lambdaMax).
3240 // State: displacements u, momenta p (kg m/s), and W, the energy (J)
3341 // dissipated so far (so E + W is conserved).
3442
@@ -46,10 +54,12 @@ export interface Params {
4654 rise: number; // m, height of the curved rest shape; 0 = straight
4755 a: number; // m, rest spacing (for drawing only; it does not enter the equations)
4856 tau: number; // s, amplitude decay time from drag (all modes); 0 = no damping
57+ tauTop: number; // s, amplitude decay time of the top of the band from the loss alone; 0 = off
4958 }
5059
5160 // Drag rate gamma (1/s).
5261 export const dragRate = (P: Params) => (P.tau > 0 ? 2 / P.tau : 0);
62+export const isDamped = (P: Params) => P.tau > 0 || P.tauTop > 0;
5363
5464 const freeRight = (P: Params) => P.ends === "cantilever";
5565 // Springs j = 0..nS-1 join nodes j - 1 and j.
@@ -57,6 +67,19 @@ export const springCount = (P: Params) => (freeRight(P) ? P.N : P.N + 1);
5767 // Joints at masses 0..nJ-1 (each needs both neighbours).
5868 const jointCount = (P: Params) => (freeRight(P) ? P.N - 1 : P.N);
5969
70+// Largest eigenvalue of L: 4 sin^2(N pi / (2 (N + 1))) between two supports;
71+// with a free end, 4 sin^2((2N - 1) pi / (2 (2N + 1))).
72+const laplacianMax = (P: Params) =>
73+ 4 * Math.sin(freeRight(P) ? ((2 * P.N - 1) * Math.PI) / (2 * (2 * P.N + 1)) : (P.N * Math.PI) / (2 * (P.N + 1))) ** 2;
74+// Loss coefficient c (kg/s); 0 when off.
75+const lossCoefficient = (P: Params) => (P.tauTop > 0 ? (2 * P.m) / (P.tauTop * laplacianMax(P)) : 0);
76+// (L x) into out.
77+function laplacian(P: Params, x: Float64Array, out: Float64Array): void {
78+ const N = P.N;
79+ for (let j = 0; j < N; j++) out[j] = 2 * x[j] - (j > 0 ? x[j - 1] : 0) - (j + 1 < N ? x[j + 1] : 0);
80+ if (freeRight(P)) out[N - 1] -= x[N - 1]; // no spring beyond the free end
81+}
82+
6083 export interface State {
6184 u: Float64Array;
6285 p: Float64Array;
@@ -165,13 +188,23 @@ export function energy(P: Params, s: State): number {
165188
166189 // Velocity Verlet (leapfrog), order 2, symplectic, as in the FPUT chain;
167190 // drag applied exactly in two half steps around it (Strang splitting),
168-// p -> p exp(-gamma h / 2), with the kinetic energy removed added to W.
191+// p -> p exp(-gamma h / 2), with the kinetic energy removed added to W. The
192+// loss, if on, likewise in two half steps, each a Crank-Nicolson step of
193+// dp/dt = -(c / m) L p (tridiagonal solve; it never adds energy, and its
194+// error is small while c q h / m << 1, true for any underdamped mode at the
195+// steps stepFor gives).
169196 export class Verlet {
170197 private P: Params;
171198 private F: Float64Array;
199+ private c: number; // loss coefficient
200+ private rhs: Float64Array; // scratch for the loss step
201+ private cp: Float64Array;
172202 constructor(P: Params, s: State) {
173203 this.P = P;
174204 this.F = new Float64Array(P.N);
205+ this.c = lossCoefficient(P);
206+ this.rhs = new Float64Array(P.N);
207+ this.cp = new Float64Array(P.N);
175208 forces(P, s.u, this.F);
176209 }
177210 // Call after changing u from outside.
@@ -188,12 +221,38 @@ export class Verlet {
188221 const { N, m } = this.P;
189222 const gamma = dragRate(this.P);
190223 if (gamma > 0) this.drag(s, Math.exp(-0.5 * gamma * h));
224+ if (this.c > 0) this.loss(s, h / 2);
191225 for (let j = 0; j < N; j++) p[j] += 0.5 * h * F[j];
192226 for (let j = 0; j < N; j++) u[j] += (h * p[j]) / m;
193227 forces(this.P, u, F);
194228 for (let j = 0; j < N; j++) p[j] += 0.5 * h * F[j];
229+ if (this.c > 0) this.loss(s, h / 2);
195230 if (gamma > 0) this.drag(s, Math.exp(-0.5 * gamma * h));
196231 }
232+ // (I + a L) p' = (I - a L) p, a = c dt / 2m, by the Thomas algorithm.
233+ private loss(s: State, dt: number): void {
234+ const { N, m } = this.P;
235+ const p = s.p, rhs = this.rhs, cp = this.cp;
236+ const a = (this.c * dt) / (2 * m);
237+ let T0 = 0;
238+ for (let j = 0; j < N; j++) T0 += p[j] * p[j];
239+ laplacian(this.P, p, rhs);
240+ for (let j = 0; j < N; j++) rhs[j] = p[j] - a * rhs[j];
241+ // Matrix I + a L: diagonal 1 + 2a (1 + a at a free end), off-diagonals -a.
242+ const diag = (j: number) => 1 + (freeRight(this.P) && j === N - 1 ? a : 2 * a);
243+ let b = diag(0);
244+ cp[0] = -a / b;
245+ p[0] = rhs[0] / b;
246+ for (let j = 1; j < N; j++) {
247+ b = diag(j) + a * cp[j - 1];
248+ cp[j] = -a / b;
249+ p[j] = (rhs[j] + a * p[j - 1]) / b;
250+ }
251+ for (let j = N - 2; j >= 0; j--) p[j] -= cp[j] * p[j + 1];
252+ let T1 = 0;
253+ for (let j = 0; j < N; j++) T1 += p[j] * p[j];
254+ s.W += (T0 - T1) / (2 * m);
255+ }
197256 private drag(s: State, f: number): void {
198257 let p2 = 0;
199258 for (let j = 0; j < s.p.length; j++) {
@@ -319,6 +378,21 @@ export function nonlinearScale(P: Params, modes: Modes): number {
319378 return V4 > 0 ? (V2 * V2) / V4 : Infinity; // V2 A^2
320379 }
321380
381+// Amplitude decay time (s) of each linear mode from drag and loss together
382+// (Infinity if neither): 1 / (1 / tau + q_n / (lambdaMax tauTop)), q_n =
383+// S_n L S_n. Exact for drag; for the loss it treats L as diagonal in the
384+// modes (exact for pinned ends, close otherwise).
385+export function decayTimes(P: Params, modes: Modes): Float64Array {
386+ const Ls = new Float64Array(P.N), lmax = laplacianMax(P);
387+ return Float64Array.from(modes.S, (row) => {
388+ laplacian(P, row, Ls);
389+ let q = 0;
390+ for (let j = 0; j < P.N; j++) q += row[j] * Ls[j];
391+ const rate = (P.tau > 0 ? 1 / P.tau : 0) + (P.tauTop > 0 ? q / (lmax * P.tauTop) : 0);
392+ return rate > 0 ? 1 / rate : Infinity;
393+ });
394+}
395+
322396 // Integrator step (s) for motion with energy E: 0.03 / omega, with omega from
323397 // the stiffest the beam can be. Linear part: the top mode, Kmax. Stretching
324398 // adds at most 12 k3 d^2 + 24 k3 D d with d = (4 E / k3)^(1/4), the largest
src/beam/index.htmlmodified+2−1View file
@@ -69,6 +69,7 @@
6969 <label>h <input id="p-rise" type="number" step="any" /> <span class="unit">m (rise of the rest shape, 0 = straight)</span></label>
7070 <label>a <input id="p-a" type="number" step="any" /> <span class="unit">m (spacing)</span></label>
7171 <label>τ <input id="p-tau" type="number" step="any" /> <span class="unit">s (decay, 0 = off)</span></label>
72+ <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>
7273 <label class="wide">stretching <select id="p-stretch">
7374 <option value="local">local (each spring)</option>
7475 <option value="kirchhoff">Kirchhoff (one shared tension)</option>
@@ -139,7 +140,7 @@
139140 </template>
140141 <template id="help-params">
141142 <h2>Parameters</h2>
142- <p>Mass, tension k, bending κ (either may be 0, not both), stretching k3, the rest spacing a (drawing only), and the decay time τ of the amplitude (drag on the masses; 0 for none). Apply restarts the beam at rest.</p>
143+ <p>Mass, tension k, bending κ (either may be 0, not both), stretching k3, the rise h of a curved rest shape (0 for straight), the rest spacing a (drawing only), the decay time τ of the amplitude (drag on the masses, all modes alike; 0 for none), and τ top, the decay time a frequency-dependent loss (dashpots alongside the springs) gives the highest modes, less for lower ones (0 for none); the derived values list each mode's decay time. Apply restarts the beam at rest.</p>
143144 <p>Stretching "local" stiffens each spring by its own stretch, as in the FPUT chain. "Kirchhoff" gives all springs one shared tension set by the total stretch, as for a beam whose ends are held apart. The two agree when every spring is stretched alike. The nonlinear scale below is the energy of mode 1 at which its stretching energy equals its linear energy.</p>
144145 </template>
145146 <template id="help-realtime">
src/beam/main.tsmodified+19−9View file
@@ -1,7 +1,7 @@
11 import "../common/style.css";
22 import "./style.css";
33 import {
4- Modes, Verlet, copyState, energy, newState, nonlinearScale, restShape, stepFor, strikeMass, strikeMode,
4+ Modes, Verlet, copyState, decayTimes, energy, isDamped, newState, nonlinearScale, restShape, stepFor, strikeMass, strikeMode,
55 type Ends, type Params, type State, type Stretch,
66 } from "./beam";
77 import { fit, indexColor } from "../common/canvas";
@@ -39,7 +39,7 @@ const pauseButton = document.getElementById("pause") as HTMLButtonElement;
3939 // clamped at both ends, with the FPUT chain's stretching coefficient.
4040 // tau: amplitude decay time from drag; 540 s makes the default strike ring
4141 // down in about 8 s of sound at the default pitch.
42-let P: Params = { N: 8, m: 1, k: 0, kappa: 20, k3: 4000, stretch: "local", ends: "clamped", rise: 0, a: 0.5, tau: 540 };
42+let P: Params = { N: 8, m: 1, k: 0, kappa: 20, k3: 4000, stretch: "local", ends: "clamped", rise: 0, a: 0.5, tau: 540, tauTop: 0 };
4343 const SAMPLE_DT = 0.01; // s between plotted samples
4444
4545 // The beam. It sits at rest until the first Strike; changing parameters puts
@@ -208,6 +208,7 @@ const paramInputs = {
208208 rise: document.getElementById("p-rise") as HTMLInputElement,
209209 a: document.getElementById("p-a") as HTMLInputElement,
210210 tau: document.getElementById("p-tau") as HTMLInputElement,
211+ tauTop: document.getElementById("p-tautop") as HTMLInputElement,
211212 };
212213 const stretchSelect = document.getElementById("p-stretch") as HTMLSelectElement;
213214 const paramStatus = document.getElementById("param-status")!;
@@ -217,6 +218,14 @@ function showParams() {
217218 endsSelect.value = P.ends;
218219 showDerived();
219220 }
221+// Amplitude decay times: one for all modes (drag alone), or per mode.
222+function decayText(): string {
223+ if (!damped()) return "none (no damping)";
224+ if (P.tauTop === 0) return `${P.tau} s (amplitude, all modes)`;
225+ const d = decayTimes(P, modes), f = (x: number) => fmtDuration(x);
226+ const list = P.N <= 8 ? Array.from(d, f).join(", ") : `${f(d[0])} … ${f(d[P.N - 1])}`;
227+ return `${list} (amplitude, by mode)`;
228+}
220229 // Values that follow from the parameters (and, for the step, the last strike).
221230 function showDerived() {
222231 const hz = Array.from(modes.omega, (w) => w / (2 * Math.PI));
@@ -226,7 +235,7 @@ function showDerived() {
226235 `normal modes ${list} Hz\n` +
227236 `ratios to mode 1 ${ratios}${P.N > 6 ? ", …" : ""}; top ${(hz[P.N - 1] / hz[0]).toFixed(1)}\n` +
228237 `nonlinear scale ${fmtEnergy(nonlinearScale(P, modes))} (mode 1: stretching = linear energy)\n` +
229- `decay time ${damped() ? `${P.tau} s (amplitude, all modes)` : "none (no damping)"}\n` +
238+ `decay time ${decayText()}\n` +
230239 `integrator step ${(1000 * H).toPrecision(3)} ms (from the strike energy; ` +
231240 `${(1000 * stepFor(P, modes, 0)).toPrecision(3)} ms when weak, ` +
232241 `${(1000 * stepFor(P, modes, 10 ** Number(energyInput.max))).toPrecision(3)} ms at the top of the slider)`;
@@ -259,14 +268,15 @@ setupMute(document.getElementById("mute") as HTMLButtonElement, player);
259268 document.getElementById("volume")!.addEventListener("input", (e) =>
260269 player.setVolume(Number((e.target as HTMLInputElement).value)),
261270 );
262-const damped = () => P.tau > 0;
271+const damped = () => isDamped(P);
263272 const pitch = () => 10 ** Number(pitchInput.value); // Hz of the lowest linear mode
264273 const speedup = () => (pitch() * 2 * Math.PI) / modes.omega[0];
265274 const modePitches = () => Array.from(modes.omega, (w) => (w / (2 * Math.PI)) * speedup());
266-// Energy decays at about 2/tau, so a drop by END_RATIO takes about
267-// (tau/2) ln(1/END_RATIO) of physical time (s of sound here).
275+// The slowest mode's energy decays at 2 / tau1 (tau1 its amplitude decay
276+// time), so a drop by END_RATIO takes about (tau1 / 2) ln(1 / END_RATIO) of
277+// physical time (s of sound here).
268278 const predictedNoteLength = () =>
269- damped() ? ((P.tau / 2) * Math.log(1 / END_RATIO)) / speedup() : UNDAMPED_LENGTH;
279+ damped() ? ((Math.max(...decayTimes(P, modes)) / 2) * Math.log(1 / END_RATIO)) / speedup() : UNDAMPED_LENGTH;
270280 function showPitch() {
271281 const f = modePitches();
272282 const top = f[f.length - 1];
@@ -461,10 +471,10 @@ endsSelect.addEventListener("change", () => setParams({ ...P, ends: endsSelect.v
461471 document.getElementById("param-apply")!.addEventListener("click", () => {
462472 const v = Object.fromEntries(Object.entries(paramInputs).map(([key, el]) => [key, Number(el.value)]));
463473 const ok = Object.entries(v).every(
464- ([key, x]) => Number.isFinite(x) && (x > 0 || (["tau", "k", "kappa", "rise"].includes(key) && x === 0)),
474+ ([key, x]) => Number.isFinite(x) && (x > 0 || (["tau", "tauTop", "k", "kappa", "rise"].includes(key) && x === 0)),
465475 );
466476 if (!ok || !(v.k > 0 || v.kappa > 0)) {
467- paramStatus.textContent = "parameters must be positive (tau, rise, k or κ may be 0, not both k and κ)";
477+ paramStatus.textContent = "parameters must be positive (τ, τ top, rise, k or κ may be 0, not both k and κ)";
468478 return;
469479 }
470480 paramStatus.textContent = "";
src/beam/render-worker.tsmodified+2−2View file
@@ -2,7 +2,7 @@
22 // (the loop, energy and Lyapunov exponent: ../common/note-render.ts), with
33 // the mode energies, N values per diagnostic sample.
44
5-import { Modes, Verlet, copyState, energy, stepFor, type Params, type State } from "./beam";
5+import { Modes, Verlet, copyState, energy, isDamped, stepFor, type Params, type State } from "./beam";
66 import { renderNote, type Note } from "../common/note-render";
77 import { shadowOf } from "./lyapunov";
88 import { signalByKey } from "./signals";
@@ -45,7 +45,7 @@ self.onmessage = (e: MessageEvent<RenderRequest>) => {
4545 for (const x of Ek) modeE.push(x);
4646 },
4747 },
48- { speedup: r.speedup, sampleRate: r.sampleRate, damped: P.tau > 0, maxLength: r.maxLength, lyapunov: r.lyapunov },
48+ { speedup: r.speedup, sampleRate: r.sampleRate, damped: isDamped(P), maxLength: r.maxLength, lyapunov: r.lyapunov },
4949 );
5050 const result: RenderResult = { id: r.id, ...note, modeEnergy: Float32Array.from(modeE.slice(0, note.energy.length * P.N)) };
5151 (self as unknown as Worker).postMessage(result, [