Bead-sheet cymbal: 2D FPUT plate with bending and dome, chaos scan, web app
8 changed files+1276−0
.gitignoreadded+3−0View file
@@ -0,0 +1,3 @@
1+__pycache__/
2+*.npz
3+*.wav
README.mdadded+53−0View file
@@ -0,0 +1,53 @@
1+# fput-plate
2+
3+The fourth chaos toy: the FPUT chain made two-dimensional, as a toy cymbal or gong. It follows `../gong-wire`, which added bending stiffness and a dome to the FPUT chain but sounded thin. The reason is the mode count: a one-dimensional wire with plate-like dispersion has fewer and fewer modes per octave as the pitch rises, so its high band holds only a handful of widely spaced partials. A two-dimensional plate has roughly a constant number of modes per unit frequency, so each octave up holds about twice as many as the one below, and that dense top is what a cymbal's wash is made of.
4+
5+## The model (`plate.py`)
6+
7+$N \times N$ unit beads on a square grid move up and down, $w_{ij}$; the boundary ring is pinned and the bending simply supported. The sheet may be built as a shallow dome $z = 16 H_d x(1-x)y(1-y)$. The energy is
8+
9+$$\sum_{\text{links}} \left[\tfrac12 T D^2 + \beta\left(z' D + \tfrac12 D^2\right)^2\right] + \tfrac12 K_b \sum_{ij} (\Delta_h w)_{ij}^2,$$
10+
11+with $D$ the height difference across a link to one of the four neighbours, $z'$ its rest slope, and $\Delta_h$ the 5-point Laplacian. The links are the FPUT-$\beta$ springs in two directions (a height difference stretches a link by about $D^2/2$), standing in for a plate's midplane stretching; bending gives plate-like dispersion; the dome adds stiffness and the three-wave coupling of curved shells. Without a dome the modes are products of sines with $\omega^2 = T\lambda + K_b\lambda^2$, $\lambda = 4\sin^2(k\pi/2(N+1)) + 4\sin^2(l\pi/2(N+1))$, and mode energies are computed in this basis (approximately, when there is a dome). Damping is friction plus dashpots between neighbours, applied explicitly each step; time stepping is leapfrog. Checks: the force equals minus the energy gradient to 2e-11, the sine products are exact modes of the flat sheet to 4e-16, and energy is conserved to 1e-5 to 1e-4.
12+
13+One simplification to keep in mind: the links pull only vertically, while a real plate's midplane can also slide sideways and partly relax the stretching. The sheet's nonlinearity is therefore stiffer than a plate's; FPUT makes the same simplification.
14+
15+## Results (2026-10-05)
16+
17+With 24 x 24 beads there are 576 modes. With bending stiffness 25 and the lowest mode at 100 Hz they reach 16.8 kHz; with bending 4 they reach only 8.6 kHz, but even then the count per octave rises steadily: 2, 7, 18, 36, 104 and 300 modes in the octaves from 100 Hz up to 6.4 kHz.
18+
19+**Exploration, undamped** (`experiments.py`, `results/explore.png`; strike centered on bead (8, 10), radius 2). Share of the energy above 4 kHz at 10, 50 and 300 ms, dome 10:
20+
21+| bending | top mode | E = 3000 | E = 10000 | E = 30000 |
22+| --- | --- | --- | --- | --- |
23+| 8 | 11.4 kHz | 5%, 25%, 60% | 30%, 59%, 69% | 50%, 70%, 74% |
24+| 25 | 16.8 kHz | 0.2%, 1%, 13% | 3%, 13%, 49% | 21%, 54%, 74% |
25+
26+The sheet needs much harder strikes than the wire or the chain (around $E = 10^4$ rather than 10 to 50), because the bending stiffness dominates the linear forces and the stretching only competes at larger amplitudes. At bending 25 and $E = 10^4$ the energy climbs into the band above 4 kHz over 300 ms, a gradual bloom; at $E = 3 \times 10^4$ it climbs within tens of milliseconds. A dome of height 10 adds a few percent; without one the picture is the same.
27+
28+**Bloom without chaos, unless damping is light.** In the web app (lowest mode ringing 4 s, highest 0.8 s) the Gong preset ($E = 10^4$) puts 22% of the energy above 4 kHz at 145 ms, and the Crash preset ($3 \times 10^4$) 47% at 152 ms. With the highest mode ringing only 0.25 s, even the Crash strike was essentially regular: a perturbation grew by about $e^{1.2}$ in all, and strikes differing by $10^{-6}$ ended identical. Strikes differing by $10^{-3}$ still drifted apart, through amplitude-dependent pitch shifts rather than exponential separation. Undamped, the same strike is chaotic (growth about $e^{12}$ in 0.6 s). So in this toy, as in the gong wire, the high-frequency bloom is robust and reproducible, while chaos needs both a hard strike and light damping of the high modes.
29+
30+**Where the chaos is** (`chaos_scan.py`). Pairs of strikes differing by $10^{-6}$ in strength, 1.5 s of sound, lowest mode at 100 Hz ringing 6 s, strike at bead (9, 10) with radius 1.6. Reported: RMS difference of the final energy share per band, in dB, for the highest mode ringing 0.8 / 2 / 5 s. Values below about 0.1 dB mean the strikes ended the same; above 1 dB they sound different.
31+
32+| bending, dome | E = 10^4 | E = 3 x 10^4 | E = 10^5 |
33+| --- | --- | --- | --- |
34+| 8, flat | 0.00 / 0.05 / 0.30 | 1.5 / 2.3 / 3.2 | 2.9 / 2.8 / 4.0 |
35+| 8, dome 10 | 0.05 / 0.63 / 1.8 | 2.2 / 4.7 / 2.3 | 3.1 / 5.5 / 2.7 |
36+| 25, flat | 0.03 / 0.02 / 0.04 | 0.02 / 0.04 / 0.04 | 0.11 / 2.5 / 2.5 |
37+| 25, dome 10 | 0.01 / 0.01 / 0.02 | 0.03 / 0.35 / 0.24 | 1.6 / 3.6 / 4.2 |
38+
39+Softer bending makes chaos far easier to reach: at bending 8 every strike of $3 \times 10^4$ or more is chaotic, while bending 25 needs about $10^5$ and a high-mode ring of 2 s or more. Presumably this is because the stretching competes with a weaker linear restoring force, so the same strike is relatively more nonlinear. A dome helps in most cases, and so does lighter damping of the high modes, with some exceptions (bending 8, dome 10, $E = 3 \times 10^4$ differs less at 5 s than at 2 s; with only one pair per setting these numbers carry noise). The app's original presets (bending 25, $E \le 3 \times 10^4$, high ring 0.8 s) all sit in the regular corner, which is why hits there sounded identical.
40+
41+**Time step.** The app steps at about 0.6 of leapfrog's linear stability limit. Waveforms differ from runs with 4 times smaller steps (the top modes accumulate phase error over hundreds of periods), but the share of energy above 4 kHz over time and the final spectrum by band agree within a few percent, which is what the sound depends on.
42+
43+## Web application
44+
45+`web/index.html` (published as a private artifact at https://claude.ai/artifact/UeCkZbQqcbzpk3e4kvCmr5) runs 8 strikes and a separate Lyapunov run in parallel Web Workers. Controls: strike energy, bending stiffness, dome height, beads per side, strike position (sliders, or click the sheet) and radius, listener position (sliders, or shift-click) or the whole sheet summed, the difference between strikes, ring times, pitch of the lowest mode (held so the top stays under 20 kHz), and length. Panels: the sheet from above with heights as color (with slow motion), strike A's energy per frequency band over time, the share above 1 and 4 kHz for all strikes, their final spectra by band, the local Lyapunov exponent with the energy, and a spectrogram of any strike. Presets: Soft ($10^3$), Gong ($10^4$), Crash ($3 \times 10^4$), Flat (no dome), all regular; and two chosen from the scan to be chaotic: Chaotic gong (bending 8, $E = 3 \times 10^4$, high ring 2 s, lowest mode 150 Hz; final spectra differ by 1.5 dB in the browser) and Chaotic crash (bending 25, $E = 10^5$, high ring 2 s; 1.1 dB). Computing a preset takes about 5 to 10 s on a 12-core laptop; the Lyapunov run, which also integrates the tangent dynamics, takes about three times as long as a strike and sets the wait. The step code works on zero-padded arrays with the substeps fused into three passes; it gives bit-for-bit the same output as the plain version and runs about 1.5 times faster. Locally: `cd web && python -m http.server`.
46+
47+## Running
48+
49+Experiments run on the office workstation (`~/work/fput-plate`, Python 3.9 with numpy), plots locally:
50+
51+ WORKERS=12 python experiments.py # results/explore.npz (weaker strikes)
52+ TAG=explore2 WORKERS=12 python experiments.py # results/explore2.npz (harder strikes)
53+ python plots.py # results/explore.png
chaos_scan.pyadded+80−0View file
@@ -0,0 +1,80 @@
1+"""Where does the struck sheet become chaotic? Pairs of strikes differing by 1e-6 in strength, with the web app's
2+damping and time stepping; how different are they after a second?
3+
4+ WORKERS=32 python chaos_scan.py # results/chaos_scan.npz and a printed table
5+
6+For each setting: bending Kb, dome Hd, strike energy E, ring time of the highest mode (lowest mode 6 s), 1.5 s of
7+sound with the lowest mode at 100 Hz. Reported: relative difference of the pickup velocity in 1.0-1.5 s, and the
8+RMS difference of the final energy share per log-spaced band (dB, bands holding > 0.1%).
9+"""
10+
11+from __future__ import annotations
12+
13+import itertools
14+import os
15+from multiprocessing import Pool
16+
17+os.environ.setdefault("OMP_NUM_THREADS", "1")
18+
19+import numpy as np # noqa: E402
20+
21+from plate import Plate, strike # noqa: E402
22+
23+HERE = os.path.dirname(os.path.abspath(__file__))
24+RES = os.path.join(HERE, "results")
25+N, SR, F1, DUR, DELTA = 24, 44100, 100.0, 1.5, 1e-6
26+
27+
28+def one(args):
29+ Kb, Hd, E, thi = args
30+ pl0 = Plate(N=N, Kb=Kb, Hd=Hd)
31+ om = pl0.omega()
32+ om1, omN = om.min(), om.max()
33+ tau = om1 / (2 * np.pi * F1)
34+ M = N + 1
35+ s1, sN = 4 * np.sin(np.pi / (2 * M)) ** 2, 4 * np.sin(N * np.pi / (2 * M)) ** 2
36+ lam1, lamN = 2 * s1, 2 * sN
37+ sLo, sHi = np.log(1000) / 6.0 * tau, np.log(1000) / thi * tau
38+ c1 = max(0.0, 2 * (sHi - sLo) / (lamN - lam1))
39+ c0 = max(0.0, 2 * sLo - c1 * lam1)
40+ pl = Plate(N=N, Kb=Kb, Hd=Hd, c0=c0, c1=c1)
41+ dtOut = 1 / SR / tau
42+ nsub = int(np.ceil(dtOut / (1.2 / omN)))
43+ h = dtOut / nsub
44+ w, v = strike(pl, E, 9, 10, 1.6, M=2, scale=[1.0, 1 + DELTA])
45+ n = int(DUR * SR)
46+ sig = np.zeros((n, 2))
47+ f = pl.force(w)
48+ for i in range(n):
49+ for _ in range(nsub):
50+ v = v + 0.5 * h * f
51+ w = w + h * v
52+ f = pl.force(w)
53+ v = v + 0.5 * h * f
54+ v = v - h * (pl.c0 * v - pl.c1 * pl.lap(v))
55+ if not np.isfinite(v).all():
56+ return args, np.nan, np.nan, True
57+ sig[i] = v[:, 15, 17]
58+ a, b = sig[int(1.0 * SR):, 0], sig[int(1.0 * SR):, 1]
59+ wave = float(np.linalg.norm(a - b) / np.linalg.norm(a))
60+ Ek = pl.mode_energies(w, v).reshape(2, -1)
61+ fk = (om / om1 * F1).ravel()
62+ edges = np.geomspace(F1, fk.max() * 1.0001, 25)
63+ bands = np.stack([Ek[:, (fk >= lo) & (fk < hi)].sum(1) for lo, hi in zip(edges[:-1], edges[1:])], 1)
64+ bands /= bands.sum(1, keepdims=True)
65+ keep = bands.mean(0) > 1e-3
66+ spec = float(np.sqrt(np.mean((10 * np.log10(bands[0, keep] / bands[1, keep])) ** 2)))
67+ return args, wave, spec, False
68+
69+
70+if __name__ == "__main__":
71+ os.makedirs(RES, exist_ok=True)
72+ grid = list(itertools.product((8.0, 25.0), (0.0, 10.0), (1e4, 3e4, 1e5), (0.8, 2.0, 5.0)))
73+ rows = []
74+ with Pool(int(os.environ.get("WORKERS", "32"))) as pool:
75+ for args, wave, spec, blew in pool.imap_unordered(one, grid):
76+ rows.append((*args, wave, spec, blew))
77+ Kb, Hd, E, thi = args
78+ print(f"Kb {Kb:4g} dome {Hd:4g} E {E:7g} high ring {thi:3g} s: "
79+ + ("BLEW UP" if blew else f"waveform diff 1-1.5 s {wave:.1e}, final band spectra differ {spec:.2f} dB"), flush=True)
80+ np.savez(os.path.join(RES, "chaos_scan.npz"), rows=np.array(rows, dtype=float))
experiments.pyadded+61−0View file
@@ -0,0 +1,61 @@
1+"""Exploration of the struck sheet (plate.py). Writes results/explore.npz; plots.py draws it.
2+
3+ python experiments.py explore # undamped: bending x dome x strike energy; energy above 1 / 4 kHz over time
4+
5+Time is in seconds of sound with the lowest mode at F_LOW Hz, as in the web app.
6+"""
7+
8+from __future__ import annotations
9+
10+import os
11+import sys
12+from multiprocessing import Pool
13+
14+os.environ.setdefault("OMP_NUM_THREADS", "1")
15+
16+import numpy as np # noqa: E402
17+
18+from plate import Plate, run, strike # noqa: E402
19+
20+HERE = os.path.dirname(os.path.abspath(__file__))
21+RES = os.path.join(HERE, "results")
22+N, CI, CJ, RAD, F_LOW = 24, 8, 10, 2.0, 100.0
23+
24+
25+def _one(args):
26+ Kb, Hd, E = args
27+ pl = Plate(N=N, Kb=Kb, Hd=Hd)
28+ om = pl.omega()
29+ tau = om.min() / (2 * np.pi * F_LOW)
30+ dt = 0.6 / np.sqrt(om.max() ** 2 + 12 * E / N)
31+ rec = max(1, int(round(0.002 / tau / dt)))
32+ w0, v0 = strike(pl, E, CI, CJ, RAD, M=2, scale=[1.0, 1 + 1e-6])
33+ r = run(pl, w0, v0, 0.3 / tau, dt, record_every=rec)
34+ f = (om / om.min() * F_LOW).ravel()
35+ Ek = r["Ek"][:, 0].reshape(len(r["t"]), -1)
36+ return Kb, Hd, E, r["t"] * tau, Ek, f, r["E"][:, 0]
37+
38+
39+def explore(pool):
40+ tag = os.environ.get("TAG", "explore")
41+ if tag == "explore2":
42+ jobs = [(Kb, Hd, E) for Kb in (8.0, 25.0) for Hd in (0.0, 10.0) for E in (3000.0, 10000.0, 30000.0)]
43+ else:
44+ jobs = [(Kb, Hd, E) for Kb in (8.0, 25.0) for Hd in (0.0, 4.0) for E in (10.0, 100.0, 1000.0)]
45+ out = {}
46+ for Kb, Hd, E, t, Ek, f, Etot in pool.imap_unordered(_one, jobs):
47+ k = f"Kb{Kb:g}_Hd{Hd:g}_E{E:g}"
48+ out[k + "_t"], out[k + "_Ek"], out[k + "_f"], out[k + "_E"] = t, Ek, f, Etot
49+ sh = Ek / Ek.sum(1, keepdims=True)
50+ h1, h4 = sh[:, f > 1000].sum(1), sh[:, f > 4000].sum(1)
51+ i10, i50 = np.argmin(abs(t - 0.01)), np.argmin(abs(t - 0.05))
52+ print(f" Kb {Kb:g} Hd {Hd:g} E {E:g}: top {f.max():.0f} Hz; above 1 kHz at 10/50/300 ms {h1[i10]:.2f} {h1[i50]:.2f} "
53+ f"{h1[-1]:.2f}; above 4 kHz {h4[i10]:.3f} {h4[i50]:.3f} {h4[-1]:.3f}; energy drift {abs(Etot[-1] / Etot[0] - 1):.1e}",
54+ flush=True)
55+ np.savez(os.path.join(RES, tag + ".npz"), **out)
56+
57+
58+if __name__ == "__main__":
59+ os.makedirs(RES, exist_ok=True)
60+ with Pool(int(os.environ.get("WORKERS", "12"))) as pool:
61+ explore(pool)
plate.pyadded+126−0View file
@@ -0,0 +1,126 @@
1+"""A struck sheet of beads: the FPUT chain made two-dimensional, as a toy cymbal or gong.
2+
3+N x N unit beads on a square grid move up and down, w_ij; the boundary ring is pinned (w = 0 on it, and the
4+bending is simply supported). The sheet may be built as a shallow dome of height Hd, rest shape
5+z = 16 Hd x (1 - x) y (1 - y). Energy:
6+ links (to the 4 neighbours, including the pinned ring): sum [ T D^2 / 2 + beta (z' D + D^2 / 2)^2 ]
7+ bending: Kb / 2 sum (Lap w)^2 (5-point Laplacian)
8+with D the height difference across a link and z' its rest slope. The links are the 2D version of the FPUT-beta
9+springs (a height difference D stretches a link by about D^2 / 2); bending gives plate-like dispersion; the dome adds
10+the quadratic (three-wave) coupling of curved shells. Without a dome the normal modes are products of sines,
11+ sin(k pi i / (N + 1)) sin(l pi j / (N + 1)), omega^2 = T lam + Kb lam^2,
12+ lam = 4 sin^2(k pi / (2 (N + 1))) + 4 sin^2(l pi / (2 (N + 1))),
13+and the mode energies below use this basis (approximate when there is a dome). Damping: friction c0 and
14+neighbour dashpots c1 (Laplacian of the velocity), applied explicitly each step. Time stepping: leapfrog.
15+"""
16+
17+from __future__ import annotations
18+
19+from dataclasses import dataclass
20+
21+import numpy as np
22+
23+
24+@dataclass(frozen=True)
25+class Plate:
26+ N: int = 24
27+ T: float = 1.0
28+ beta: float = 1.0
29+ Kb: float = 4.0
30+ Hd: float = 0.0
31+ c0: float = 0.0
32+ c1: float = 0.0
33+
34+ def _z(self):
35+ g = np.arange(self.N + 2) / (self.N + 1)
36+ return 16 * self.Hd * np.outer(g * (1 - g), g * (1 - g))
37+
38+ @property
39+ def zh(self):
40+ """Rest slopes of the horizontal links, (N, N + 1)."""
41+ return np.diff(self._z()[1:-1, :], axis=1)
42+
43+ @property
44+ def zv(self):
45+ """Rest slopes of the vertical links, (N + 1, N)."""
46+ return np.diff(self._z()[:, 1:-1], axis=0)
47+
48+ def lam(self):
49+ s = 4 * np.sin(np.arange(1, self.N + 1) * np.pi / (2 * (self.N + 1))) ** 2
50+ return s[:, None] + s[None, :]
51+
52+ def omega(self):
53+ lam = self.lam()
54+ return np.sqrt(self.T * lam + self.Kb * lam**2)
55+
56+ def sine(self):
57+ N = self.N
58+ j = np.arange(1, N + 1)
59+ return np.sqrt(2 / (N + 1)) * np.sin(np.outer(j, j) * np.pi / (N + 1))
60+
61+ @staticmethod
62+ def lap(u):
63+ """5-point Laplacian of (M, N, N) with zero boundary."""
64+ p = np.pad(u, ((0, 0), (1, 1), (1, 1)))
65+ return p[:, 2:, 1:-1] + p[:, :-2, 1:-1] + p[:, 1:-1, 2:] + p[:, 1:-1, :-2] - 4 * u
66+
67+ def force(self, w):
68+ p = np.pad(w, ((0, 0), (1, 1), (1, 1)))
69+ Dh = np.diff(p[:, 1:-1, :], axis=2) # (M, N, N + 1)
70+ Dv = np.diff(p[:, :, 1:-1], axis=1) # (M, N + 1, N)
71+ zh, zv, T, b = self.zh, self.zv, self.T, self.beta
72+ eh = zh * Dh + 0.5 * Dh * Dh
73+ ev = zv * Dv + 0.5 * Dv * Dv
74+ fh = T * Dh + 2 * b * eh * (zh + Dh)
75+ fv = T * Dv + 2 * b * ev * (zv + Dv)
76+ f = (fh[:, :, 1:] - fh[:, :, :-1]) + (fv[:, 1:, :] - fv[:, :-1, :])
77+ if self.Kb:
78+ f -= self.Kb * self.lap(self.lap(w))
79+ return f
80+
81+ def energy(self, w, v):
82+ p = np.pad(w, ((0, 0), (1, 1), (1, 1)))
83+ Dh = np.diff(p[:, 1:-1, :], axis=2)
84+ Dv = np.diff(p[:, :, 1:-1], axis=1)
85+ eh = self.zh * Dh + 0.5 * Dh * Dh
86+ ev = self.zv * Dv + 0.5 * Dv * Dv
87+ U = (0.5 * self.T * Dh**2 + self.beta * eh**2).sum((1, 2)) + (0.5 * self.T * Dv**2 + self.beta * ev**2).sum((1, 2))
88+ U += 0.5 * self.Kb * (self.lap(w) ** 2).sum((1, 2))
89+ return 0.5 * (v * v).sum((1, 2)) + U
90+
91+ def mode_energies(self, w, v):
92+ S = self.sine()
93+ q = np.einsum("ai,mij,jb->mab", S, w, S)
94+ p = np.einsum("ai,mij,jb->mab", S, v, S)
95+ return 0.5 * (p * p + (self.omega() * q) ** 2)
96+
97+
98+def strike(pl: Plate, E, ci, cj, radius, M=1, scale=None):
99+ i = np.arange(1, pl.N + 1)
100+ prof = np.exp(-0.5 * (((i[:, None] - ci) ** 2 + (i[None, :] - cj) ** 2) / radius**2))
101+ prof *= np.sqrt(2 * E / (prof**2).sum())
102+ v = np.tile(prof, (M, 1, 1))
103+ if scale is not None:
104+ v *= np.asarray(scale)[:, None, None]
105+ return np.zeros_like(v), v
106+
107+
108+def run(pl: Plate, w0, v0, t_end, dt, record_every=100):
109+ w, v = np.array(w0, float), np.array(v0, float)
110+ n = int(round(t_end / dt))
111+ rec = dict(t=[], Ek=[], E=[])
112+ f = pl.force(w)
113+ for i in range(n):
114+ v = v + 0.5 * dt * f
115+ w = w + dt * v
116+ f = pl.force(w)
117+ v = v + 0.5 * dt * f
118+ if pl.c0 or pl.c1:
119+ v = v - dt * (pl.c0 * v - pl.c1 * pl.lap(v))
120+ if (i + 1) % record_every == 0:
121+ rec["t"].append((i + 1) * dt)
122+ rec["Ek"].append(pl.mode_energies(w, v))
123+ rec["E"].append(pl.energy(w, v))
124+ out = {k: np.array(x) for k, x in rec.items()}
125+ out.update(w=w, v=v)
126+ return out
plots.pyadded+41−0View file
@@ -0,0 +1,41 @@
1+"""Figure from results/explore2.npz (written by experiments.py with TAG=explore2). Run: python plots.py"""
2+
3+from __future__ import annotations
4+
5+import os
6+
7+import matplotlib
8+
9+matplotlib.use("Agg")
10+import matplotlib.pyplot as plt # noqa: E402
11+import numpy as np # noqa: E402
12+
13+HERE = os.path.dirname(os.path.abspath(__file__))
14+RES = os.path.join(HERE, "results")
15+
16+if __name__ == "__main__":
17+ d = np.load(os.path.join(RES, "explore2.npz"))
18+ fig, ax = plt.subplots(2, 2, figsize=(14, 8))
19+ for j, Kb in enumerate((8.0, 25.0)):
20+ k = f"Kb{Kb:g}_Hd10_E{10000:g}"
21+ t, Ek, f = d[k + "_t"], d[k + "_Ek"], d[k + "_f"]
22+ order = np.argsort(f)
23+ sh = Ek / Ek.sum(1, keepdims=True)
24+ edges = np.geomspace(f.min(), f.max() * 1.0001, 37)
25+ B = np.stack([sh[:, (f >= a) & (f < b)].sum(1) for a, b in zip(edges[:-1], edges[1:])], 1)
26+ ax[0, j].pcolormesh(t * 1000, np.sqrt(edges[:-1] * edges[1:]), np.log10(B.T + 1e-12), vmin=-6, vmax=0, cmap="magma",
27+ shading="auto", rasterized=True)
28+ ax[0, j].set(yscale="log", title=f"bending {Kb:g}, dome 10, E = 10000: energy share per band", xlabel="time (ms of sound)",
29+ ylabel="pitch (Hz)")
30+ for E, col in zip((3000.0, 10000.0, 30000.0), ("C0", "C1", "C3")):
31+ k = f"Kb{Kb:g}_Hd10_E{E:g}"
32+ t, Ek, f = d[k + "_t"], d[k + "_Ek"], d[k + "_f"]
33+ sh = Ek / Ek.sum(1, keepdims=True)
34+ ax[1, j].plot(t * 1000, sh[:, f > 4000].sum(1), color=col, label=f"E = {E:g}, above 4 kHz")
35+ ax[1, j].plot(t * 1000, sh[:, f > 1000].sum(1), color=col, lw=0.8, ls="--", label=f"E = {E:g}, above 1 kHz")
36+ ax[1, j].set(ylim=(0, 1), xlabel="time (ms of sound)", ylabel="share of energy", title=f"bending {Kb:g}, dome 10, undamped")
37+ ax[1, j].grid(alpha=0.3)
38+ ax[1, 0].legend(fontsize=7)
39+ fig.tight_layout()
40+ fig.savefig(os.path.join(RES, "explore.png"), dpi=90)
41+ print("wrote results/explore.png")
results/explore.pngadded+0−0View file
Binary file not shown.
web/index.htmladded+912−0View file
@@ -0,0 +1,912 @@
1+<meta charset="utf-8">
2+<title>Bead Sheet Cymbal</title>
3+<link rel="preconnect" href="https://fonts.googleapis.com">
4+<link rel="preconnect" href="https://fonts.gstatic.com" crossorigin>
5+<link rel="stylesheet" href="https://fonts.googleapis.com/css2?family=Archivo:wdth,wght@75..100,400..800&family=IBM+Plex+Mono:wght@400;500&display=swap">
6+<style>
7+:root {
8+ --paper: #eef2f5;
9+ --panel: #f8fafb;
10+ --ink: #172334;
11+ --muted: #5b6b7d;
12+ --rule: #cbd5df;
13+ --grid: #dde4eb;
14+ --x: #2563a6;
15+ --y: #a23b72;
16+ --chaos: #e8b931;
17+ --brass: #b07a22;
18+ --up: #c2410c;
19+ --down: #2563a6;
20+ --chaos-wash: rgba(232, 185, 49, 0.18);
21+ --faint: rgba(23, 35, 52, 0.16);
22+ --focus: #2563a6;
23+ --sans: "Archivo", "Helvetica Neue", Arial, sans-serif;
24+ --mono: "IBM Plex Mono", ui-monospace, "SFMono-Regular", Menlo, monospace;
25+}
26+@media (prefers-color-scheme: dark) {
27+ :root:not([data-theme="light"]) {
28+ color-scheme: dark;
29+ --paper: #0f161f; --panel: #16202c; --ink: #e3e9f0; --muted: #93a3b5; --rule: #2a3747; --grid: #1f2a37;
30+ --x: #6fa8e8; --y: #d57bae; --chaos: #e8b931; --brass: #d8a64f; --up: #f08a5d; --down: #6fa8e8;
31+ --chaos-wash: rgba(232, 185, 49, 0.14); --faint: rgba(227, 233, 240, 0.16); --focus: #6fa8e8;
32+ }
33+}
34+:root[data-theme="dark"] {
35+ color-scheme: dark;
36+ --paper: #0f161f; --panel: #16202c; --ink: #e3e9f0; --muted: #93a3b5; --rule: #2a3747; --grid: #1f2a37;
37+ --x: #6fa8e8; --y: #d57bae; --chaos: #e8b931; --brass: #d8a64f; --up: #f08a5d; --down: #6fa8e8;
38+ --chaos-wash: rgba(232, 185, 49, 0.14); --faint: rgba(227, 233, 240, 0.16); --focus: #6fa8e8;
39+}
40+* { box-sizing: border-box; }
41+body { background: var(--paper); color: var(--ink); font-family: var(--sans); font-size: 15px; line-height: 1.5; padding-inline: 16px; padding-block: 20px 48px; }
42+.wrap { max-width: 1280px; margin: 0 auto; display: grid; gap: 20px; }
43+header { display: grid; gap: 6px; }
44+h1 { font-size: clamp(26px, 4vw, 38px); font-weight: 750; font-stretch: 78%; letter-spacing: -0.01em; line-height: 1.05; margin: 0; text-wrap: balance; }
45+h1 .eq { font-family: var(--mono); font-weight: 400; font-size: 0.42em; color: var(--muted); letter-spacing: 0; display: block; margin-top: 8px; font-stretch: 100%; }
46+.lede { max-width: 78ch; color: var(--muted); margin: 0; }
47+.layout { display: grid; grid-template-columns: 300px minmax(0, 1fr); gap: 20px; align-items: start; }
48+@media (max-width: 900px) { .layout { grid-template-columns: 1fr; } }
49+aside { background: var(--panel); border: 1px solid var(--rule); border-radius: 6px; padding: 16px; display: grid; grid-template-columns: minmax(0, 1fr); gap: 16px; position: sticky; top: calc(env(safe-area-inset-top, 0px) + 12px); }
50+@media (max-width: 900px) { aside { position: static; } }
51+.label { font-size: 11px; font-weight: 650; letter-spacing: 0.08em; text-transform: uppercase; color: var(--muted); }
52+.presets { display: flex; flex-wrap: wrap; gap: 6px; }
53+.chip { font: inherit; font-size: 13px; padding: 4px 10px; border-radius: 999px; border: 1px solid var(--rule); background: var(--paper); color: var(--ink); cursor: pointer; }
54+.chip:hover { border-color: var(--muted); }
55+.chip[aria-pressed="true"] { background: var(--ink); color: var(--paper); border-color: var(--ink); }
56+.field { display: grid; gap: 4px; }
57+.field .row { display: flex; justify-content: space-between; align-items: baseline; gap: 8px; }
58+.field label { font-size: 14px; font-weight: 550; }
59+.val { font-family: var(--mono); font-size: 13px; font-variant-numeric: tabular-nums; color: var(--ink); white-space: nowrap; }
60+.hint { font-size: 12px; color: var(--muted); margin: 0; }
61+input[type="range"] { width: 100%; accent-color: var(--x); }
62+input[type="range"]:disabled { opacity: 0.4; }
63+.two { display: grid; grid-template-columns: 1fr 1fr; gap: 8px; }
64+.player { display: grid; gap: 8px; }
65+.btns { display: grid; grid-template-columns: 1fr 1fr; gap: 8px; }
66+.btn { font: inherit; font-weight: 600; font-size: 14px; padding: 9px 10px; border-radius: 6px; border: 1px solid var(--rule); background: var(--paper); color: var(--ink); cursor: pointer; display: flex; align-items: center; justify-content: center; gap: 8px; }
67+.btn:hover { border-color: var(--muted); }
68+.btn .sw { width: 10px; height: 10px; border-radius: 2px; flex: none; }
69+.btn[aria-pressed="true"] { background: var(--ink); color: var(--paper); border-color: var(--ink); }
70+button:focus-visible, input:focus-visible, canvas:focus-visible { outline: 2px solid var(--focus); outline-offset: 2px; }
71+.stats { display: grid; grid-template-columns: auto 1fr; gap: 4px 12px; font-size: 13px; margin: 0; }
72+.stats dt { color: var(--muted); }
73+.stats dd { margin: 0; font-family: var(--mono); font-variant-numeric: tabular-nums; text-align: right; }
74+.status { font-size: 12px; color: var(--muted); min-height: 1.5em; font-family: var(--mono); }
75+main { display: grid; gap: 16px; min-width: 0; }
76+.pair { display: grid; grid-template-columns: 1fr 1fr; gap: 16px; }
77+@media (max-width: 640px) { .pair { grid-template-columns: 1fr; } }
78+figure { margin: 0; display: grid; gap: 6px; min-width: 0; }
79+figcaption { display: grid; gap: 2px; }
80+figcaption b { font-weight: 650; font-size: 15px; }
81+figcaption span { font-size: 13px; color: var(--muted); max-width: 82ch; }
82+.plot { position: relative; width: 100%; background: var(--panel); border: 1px solid var(--rule); border-radius: 6px; }
83+.plot canvas { position: absolute; inset: 0; width: 100%; height: 100%; display: block; }
84+.sheetwrap { display: grid; grid-template-columns: minmax(0, 360px) minmax(0, 1fr); gap: 16px; align-items: start; }
85+@media (max-width: 760px) { .sheetwrap { grid-template-columns: 1fr; } }
86+.sheet { aspect-ratio: 1 / 1; max-width: 100%; cursor: crosshair; }
87+.wide { aspect-ratio: 3.2 / 1; min-height: 170px; }
88+.half { aspect-ratio: 1.6 / 1; min-height: 170px; }
89+.short { aspect-ratio: 4.6 / 1; min-height: 140px; }
90+.setuptext { display: grid; gap: 10px; font-size: 14px; }
91+.setuptext p { margin: 0; max-width: 70ch; }
92+.setuptext .btn { justify-self: start; }
93+.strikes { display: flex; flex-wrap: wrap; gap: 6px; align-items: center; }
94+.sk { font: inherit; font-family: var(--mono); font-size: 13px; min-width: 38px; padding: 4px 8px; border-radius: 6px; border: 1px solid var(--rule); background: var(--paper); color: var(--ink); cursor: pointer; display: inline-flex; align-items: center; gap: 6px; justify-content: center; }
95+.sk .sw { width: 8px; height: 8px; border-radius: 2px; }
96+.sk[aria-pressed="true"] { background: var(--ink); color: var(--paper); border-color: var(--ink); }
97+.sk:disabled { opacity: 0.4; cursor: default; }
98+.strikes .btn { padding: 4px 12px; font-size: 13px; margin-left: auto; }
99+.legend { display: flex; flex-wrap: wrap; gap: 14px; font-size: 12px; color: var(--muted); }
100+.legend i { display: inline-block; width: 18px; height: 3px; vertical-align: middle; margin-right: 6px; border-radius: 2px; }
101+footer { color: var(--muted); font-size: 13px; max-width: 82ch; }
102+footer p { margin: 0 0 8px; }
103+/* keep the control column inside its panel whatever the readouts say */
104+aside > * { min-width: 0; }
105+aside .field .row { flex-wrap: wrap; }
106+aside .val { white-space: normal; text-align: right; overflow-wrap: anywhere; }
107+</style>
108+
109+<div class="wrap">
110+ <header>
111+ <h1>Bead Sheet Cymbal
112+ <span class="eq">a square sheet of beads: stretching links + bending + dome, pinned edges, struck once</span>
113+ </h1>
114+ <p class="lede">The FPUT chain made two-dimensional. Beads on a square sheet move up and down; the links between neighbours stiffen when stretched (which couples the modes), the sheet resists bending (which gives it a plate's spectrum, about 550 modes between 100 Hz and 17 kHz, denser toward the top), and it can be built as a dome. Unlike a wire, a sheet has more and more modes per octave as the pitch rises, so a hard strike can spread its energy into a dense high-frequency wash, as a cymbal or gong does.</p>
115+ </header>
116+
117+ <div class="layout">
118+ <aside aria-label="Controls">
119+ <div class="player">
120+ <div class="label">Listen</div>
121+ <div class="btns">
122+ <button class="btn" id="playA" type="button"><span class="sw" style="background: var(--x)"></span>Strike A</button>
123+ <button class="btn" id="playB" type="button"><span class="sw" style="background: var(--y)"></span>Strike B</button>
124+ <button class="btn" id="playAB" type="button">A, then B</button>
125+ <button class="btn" id="stop" type="button">Stop</button>
126+ </div>
127+ <div class="status" id="status">Computing…</div>
128+ </div>
129+
130+ <div class="field">
131+ <div class="row"><label>Listen at</label><span class="val" id="lisv"></span></div>
132+ <div class="presets" role="group" aria-label="Listener">
133+ <button class="chip" type="button" id="lisOne" aria-pressed="true">One point</button>
134+ <button class="chip" type="button" id="lisSum" aria-pressed="false">Whole sheet (sum)</button>
135+ </div>
136+ <div class="two">
137+ <input id="lx" type="range" min="0.05" max="0.95" step="0.01" aria-label="Listener across">
138+ <input id="ly" type="range" min="0.05" max="0.95" step="0.01" aria-label="Listener down">
139+ </div>
140+ <p class="hint">Shift-click the sheet to move the listener. The sum is roughly what a distant listener hears.</p>
141+ </div>
142+
143+ <div class="field">
144+ <div class="label">Presets</div>
145+ <div class="presets">
146+ <button class="chip" type="button" data-preset="soft">Soft</button>
147+ <button class="chip" type="button" data-preset="gong">Gong</button>
148+ <button class="chip" type="button" data-preset="crash">Crash</button>
149+ <button class="chip" type="button" data-preset="flat">Flat, no dome</button>
150+ <button class="chip" type="button" data-preset="cgong">Chaotic gong</button>
151+ <button class="chip" type="button" data-preset="ccrash">Chaotic crash</button>
152+ </div>
153+ </div>
154+
155+ <div class="field">
156+ <div class="row"><label for="E0">Strike energy E</label><span class="val" id="E0v"></span></div>
157+ <input id="E0" type="range" min="0" max="1000" step="1">
158+ <p class="hint">For hits that differ from one another (chaos), the strike must be hard and the high modes must ring long: with bending 8, E ≈ 30k and high modes ringing 2 s or more; with bending 25, E ≈ 100k. A dome helps. See the two Chaotic presets.</p>
159+ </div>
160+ <div class="field">
161+ <div class="row"><label for="Kb">Bending stiffness</label><span class="val" id="Kbv"></span></div>
162+ <input id="Kb" type="range" min="0" max="40" step="0.5">
163+ <p class="hint">More bending spreads the modes higher (top shown under Pitch) but takes a harder strike to set the nonlinearity going.</p>
164+ </div>
165+ <div class="field">
166+ <div class="row"><label for="Hd">Dome height</label><span class="val" id="Hdv"></span></div>
167+ <input id="Hd" type="range" min="0" max="30" step="0.5">
168+ </div>
169+ <div class="field">
170+ <div class="row"><label for="N">Beads per side</label><span class="val" id="Nv"></span></div>
171+ <input id="N" type="range" min="12" max="32" step="1">
172+ </div>
173+ <div class="field">
174+ <div class="row"><label>Strike position</label><span class="val" id="sposv"></span></div>
175+ <div class="two">
176+ <input id="sx" type="range" min="0.08" max="0.92" step="0.01" aria-label="Strike across">
177+ <input id="sy" type="range" min="0.08" max="0.92" step="0.01" aria-label="Strike down">
178+ </div>
179+ <p class="hint">Or click the sheet.</p>
180+ </div>
181+ <div class="field">
182+ <div class="row"><label for="rad">Strike radius</label><span class="val" id="radv"></span></div>
183+ <input id="rad" type="range" min="0.6" max="5" step="0.1">
184+ </div>
185+ <div class="field">
186+ <div class="row"><label for="dl">Difference between strikes</label><span class="val" id="dlv"></span></div>
187+ <input id="dl" type="range" min="-14" max="-2" step="0.5">
188+ <p class="hint">1e-6 shows chaos only. Real hits differ by about 1e-2, which makes them sound different even without chaos, through pitch shifts that depend on the amplitude.</p>
189+ </div>
190+ <div class="field">
191+ <div class="row"><label for="tlo">Ring time, lowest mode</label><span class="val" id="tlov"></span></div>
192+ <input id="tlo" type="range" min="0.5" max="10" step="0.1">
193+ </div>
194+ <div class="field">
195+ <div class="row"><label for="thi">Ring time, highest mode</label><span class="val" id="thiv"></span></div>
196+ <input id="thi" type="range" min="0.05" max="5" step="0.05">
197+ </div>
198+ <div class="field">
199+ <div class="row"><label for="f1">Pitch of the lowest mode</label><span class="val" id="f1v"></span></div>
200+ <input id="f1" type="range" min="40" max="300" step="1">
201+ <p class="hint" id="f1note" hidden>Held lower than the slider so the highest mode stays under 20 kHz.</p>
202+ </div>
203+ <div class="field">
204+ <div class="row"><label for="dur">Length</label><span class="val" id="durv"></span></div>
205+ <input id="dur" type="range" min="1" max="4" step="0.5">
206+ </div>
207+
208+ <dl class="stats" aria-live="polite">
209+ <dt>Energy above 4 kHz, peak</dt><dd id="sBloom">–</dd>
210+ <dt>… reached at</dt><dd id="sBloomT">–</dd>
211+ <dt>Chaotic phase ends</dt><dd id="sChaos">–</dd>
212+ <dt>Peak separation rate</dt><dd id="sLam">–</dd>
213+ <dt>Final spectra, 8 strikes</dt><dd id="sSpread">–</dd>
214+ </dl>
215+ </aside>
216+
217+ <main>
218+ <figure>
219+ <figcaption><b>Physical setup</b><span>The sheet seen from above. Each square is a bead; its color is its height (warm up, cool down), drawn magnified. The edges are pinned. Click to move the strike (gold); shift-click to move the listener (ring).</span></figcaption>
220+ <div class="sheetwrap">
221+ <div class="plot sheet"><canvas id="cSheet" tabindex="0" aria-label="The sheet; click to move the strike, shift-click to move the listener"></canvas></div>
222+ <div class="setuptext">
223+ <p>Each bead is tied to its four neighbours by links that pull harder as they stretch: a height difference D between neighbours stretches the link by about D²/2, so its stiffness grows with D². This is the FPUT chain's nonlinearity, now in two directions, and it plays the role of a cymbal's midplane stretching, passing energy from mode to mode.</p>
224+ <p>The sheet also resists bending, so short ripples vibrate much faster than long ones, as in a metal plate. Built as a dome, it stiffens and gains the three-wave coupling of curved shells. The strike gives a patch of beads a sudden upward velocity; damping is friction plus dashpots between neighbours, strongest for the shortest ripples.</p>
225+ <button class="btn" id="slow" type="button" aria-pressed="false">Slow motion</button>
226+ </div>
227+ </div>
228+ </figure>
229+
230+ <figure>
231+ <figcaption><b>Where the energy is</b><span>Share of strike A's energy in each frequency band (modes grouped by pitch, log scale) over time. A cascade shows as energy climbing toward the top.</span></figcaption>
232+ <div class="plot wide"><canvas id="cModes"></canvas><canvas id="cModesO"></canvas></div>
233+ </figure>
234+
235+ <div class="pair">
236+ <figure>
237+ <figcaption><b>How high it reaches</b><span>Share of each strike's energy above 4 kHz (solid) and above 1 kHz (thin), for all 8 strikes. Shaded: the chaotic phase of strike A.</span></figcaption>
238+ <div class="plot half"><canvas id="cHi"></canvas><canvas id="cHiO"></canvas></div>
239+ </figure>
240+ <figure>
241+ <figcaption><b>The spectrum it rings with</b><span>Each strike's final share of energy per frequency band.</span></figcaption>
242+ <div class="plot half"><canvas id="cFinal"></canvas></div>
243+ </figure>
244+ </div>
245+ <div class="legend"><span><i style="background: var(--x)"></i>Strike A</span><span><i style="background: var(--y)"></i>Strike B</span><span><i style="background: var(--faint)"></i>6 more strikes</span></div>
246+
247+ <figure>
248+ <figcaption><b>How fast nearby strikes separate</b><span>Local Lyapunov exponent of strike A (per second, 50 ms windows). The dashed line is the energy in dB.</span></figcaption>
249+ <div class="plot short"><canvas id="cLam"></canvas><canvas id="cLamO"></canvas></div>
250+ </figure>
251+
252+ <figure>
253+ <figcaption><b>What you hear</b><span>Spectrogram of the selected strike at the listener. Click a strike, or use the ← and → keys, to flip between them; the color scale is shared, so differences in level are real.</span></figcaption>
254+ <div class="strikes" role="group" aria-label="Choose a strike" id="strikeSel"></div>
255+ <div class="plot wide"><canvas id="cSpec"></canvas><canvas id="cSpecO"></canvas></div>
256+ </figure>
257+
258+ <footer>
259+ <p>Computed in your browser, the 8 strikes in parallel: leapfrog time stepping with the damping applied each step. Mode energies use the flat sheet's modes (products of sines), so with a dome they are approximate. Time is scaled so that the lowest mode sounds at the chosen pitch, held low enough that the highest stays below 20 kHz. The links pull only along the vertical; a real plate's midplane can also slide sideways and partly relax the stretching, so this sheet's nonlinearity is stiffer than a plate's, the same simplification FPUT makes.</p>
260+ </footer>
261+ </main>
262+ </div>
263+</div>
264+
265+<script id="worker-src" type="text/plain">
266+function setupSheet(p) {
267+ const N = p.N, M = N + 1;
268+ const s = new Float64Array(N);
269+ for (let k = 0; k < N; k++) s[k] = 4 * Math.sin((k + 1) * Math.PI / (2 * M)) ** 2;
270+ const om = new Float64Array(N * N);
271+ for (let a = 0; a < N; a++) for (let b = 0; b < N; b++) { const lam = s[a] + s[b]; om[a * N + b] = Math.sqrt(p.T * lam + p.Kb * lam * lam); }
272+ const S = new Float64Array(N * N), c = Math.sqrt(2 / M);
273+ for (let i = 0; i < N; i++) for (let k = 0; k < N; k++) S[i * N + k] = c * Math.sin((i + 1) * (k + 1) * Math.PI / M);
274+ // rest slopes of the links: z = 16 Hd x (1 - x) y (1 - y) on the padded grid
275+ const z = (i, j) => 16 * p.Hd * (i / M) * (1 - i / M) * (j / M) * (1 - j / M);
276+ const zh = new Float64Array(N * M), zv = new Float64Array(M * N);
277+ for (let i = 0; i < N; i++) for (let j = 0; j < M; j++) zh[i * M + j] = z(i + 1, j + 1) - z(i + 1, j);
278+ for (let i = 0; i < M; i++) for (let j = 0; j < N; j++) zv[i * N + j] = z(i + 1, j + 1) - z(i, j + 1);
279+ return { om, S, zh, zv };
280+}
281+
282+function strikeRun(p, sh, scale, opts) {
283+ // Arrays are padded with a ring of zeros (the pinned boundary), so the stencils need no edge tests: node (i, j)
284+ // sits at (i + 1) * P + j + 1. The steps are fused into three passes per substep; the arithmetic, and so the
285+ // result, is the same as the plain kick, drift, force, kick, damp sequence.
286+ const N = p.N, M = N + 1, P = N + 2, NN = N * N, PP = P * P, T = p.T, beta = p.beta, Kb = p.Kb, zh = sh.zh, zv = sh.zv;
287+ const c0 = p.c0, c1 = p.c1, tan = !!opts.tangent;
288+ let w = new Float64Array(PP), v = new Float64Array(PP), vb = new Float64Array(PP);
289+ const f = new Float64Array(PP), L1 = new Float64Array(PP);
290+ let dw = null, dv = null, dvb = null, g = null, dL1 = null, ch = null, cv = null;
291+ const fvA = new Float64Array(N), fvB = new Float64Array(N), gvA = new Float64Array(N), gvB = new Float64Array(N);
292+ let n2 = 0;
293+ for (let i = 0; i < N; i++) for (let j = 0; j < N; j++) {
294+ const q = Math.exp(-0.5 * ((i + 1 - p.ci) ** 2 + (j + 1 - p.cj) ** 2) / (p.rad * p.rad));
295+ v[(i + 1) * P + j + 1] = q; n2 += q * q;
296+ }
297+ const amp = Math.sqrt(2 * p.E0 / n2) * scale;
298+ for (let k = 0; k < PP; k++) v[k] *= amp;
299+ if (tan) {
300+ dw = new Float64Array(PP); dv = new Float64Array(PP); dvb = new Float64Array(PP); g = new Float64Array(PP); dL1 = new Float64Array(PP);
301+ ch = new Float64Array(N * M); cv = new Float64Array(M * N);
302+ for (let i = 0; i < N; i++) for (let j = 0; j < N; j++) { const k = i * N + j, q = (i + 1) * P + j + 1; dw[q] = Math.sin(1.7 * k + 0.3); dv[q] = Math.cos(2.3 * k + 0.1); }
303+ normalize(dw, dv);
304+ }
305+
306+ function lapInto(u, out) {
307+ for (let i = 1; i <= N; i++) { const r = i * P; for (let q = r + 1; q <= r + N; q++) out[q] = u[q - P] + u[q + P] + u[q - 1] + u[q + 1] - 4 * u[q]; }
308+ }
309+ // force on u, written to out, and the half kick vv += hh * out
310+ function forceKick(u, out, vv, hh) {
311+ for (let j = 0; j < N; j++) {
312+ const d = u[P + j + 1], zz = zv[j], e = zz * d + 0.5 * d * d;
313+ fvA[j] = T * d + 2 * beta * e * (zz + d);
314+ }
315+ lapInto(u, L1);
316+ for (let i = 0; i < N; i++) {
317+ const r = (i + 1) * P, rh = i * M, rv = (i + 1) * N;
318+ for (let j = 0; j < N; j++) {
319+ const d = u[r + P + j + 1] - u[r + j + 1], zz = zv[rv + j], e = zz * d + 0.5 * d * d;
320+ fvB[j] = T * d + 2 * beta * e * (zz + d);
321+ }
322+ let d = u[r + 1], zz = zh[rh], e = zz * d + 0.5 * d * d, hp = T * d + 2 * beta * e * (zz + d);
323+ for (let j = 0; j < N; j++) {
324+ const q = r + j + 1;
325+ d = u[q + 1] - u[q]; zz = zh[rh + j + 1]; e = zz * d + 0.5 * d * d;
326+ const hn = T * d + 2 * beta * e * (zz + d);
327+ const fo = hn - hp + fvB[j] - fvA[j] - Kb * (L1[q - P] + L1[q + P] + L1[q - 1] + L1[q + 1] - 4 * L1[q]);
328+ out[q] = fo; vv[q] += hh * fo; hp = hn;
329+ }
330+ fvA.set(fvB);
331+ }
332+ }
333+ // the same, also keeping the link stiffnesses for tforceKick
334+ function forceKickKeep(u, out, vv, hh) {
335+ for (let j = 0; j < N; j++) {
336+ const d = u[P + j + 1], zz = zv[j], e = zz * d + 0.5 * d * d;
337+ fvA[j] = T * d + 2 * beta * e * (zz + d);
338+ { const s = zz + d; cv[j] = T + 2 * beta * (s * s + zz * d + 0.5 * d * d); }
339+ }
340+ lapInto(u, L1);
341+ for (let i = 0; i < N; i++) {
342+ const r = (i + 1) * P, rh = i * M, rv = (i + 1) * N;
343+ for (let j = 0; j < N; j++) {
344+ const d = u[r + P + j + 1] - u[r + j + 1], zz = zv[rv + j], e = zz * d + 0.5 * d * d;
345+ fvB[j] = T * d + 2 * beta * e * (zz + d);
346+ { const s = zz + d; cv[rv + j] = T + 2 * beta * (s * s + zz * d + 0.5 * d * d); }
347+ }
348+ let d = u[r + 1], zz = zh[rh], e = zz * d + 0.5 * d * d, hp = T * d + 2 * beta * e * (zz + d);
349+ { const s = zz + d; ch[rh] = T + 2 * beta * (s * s + zz * d + 0.5 * d * d); }
350+ for (let j = 0; j < N; j++) {
351+ const q = r + j + 1;
352+ d = u[q + 1] - u[q]; zz = zh[rh + j + 1]; e = zz * d + 0.5 * d * d;
353+ const hn = T * d + 2 * beta * e * (zz + d);
354+ { const s = zz + d; ch[rh + j + 1] = T + 2 * beta * (s * s + zz * d + 0.5 * d * d); }
355+ const fo = hn - hp + fvB[j] - fvA[j] - Kb * (L1[q - P] + L1[q + P] + L1[q - 1] + L1[q + 1] - 4 * L1[q]);
356+ out[q] = fo; vv[q] += hh * fo; hp = hn;
357+ }
358+ fvA.set(fvB);
359+ }
360+ }
361+ // linearized force on du at the state of the last forceKickKeep, and the half kick
362+ function tforceKick(du, out, vv, hh) {
363+ for (let j = 0; j < N; j++) gvA[j] = cv[j] * du[P + j + 1];
364+ lapInto(du, dL1);
365+ for (let i = 0; i < N; i++) {
366+ const r = (i + 1) * P, rh = i * M, rv = (i + 1) * N;
367+ for (let j = 0; j < N; j++) gvB[j] = cv[rv + j] * (du[r + P + j + 1] - du[r + j + 1]);
368+ let hp = ch[rh] * du[r + 1];
369+ for (let j = 0; j < N; j++) {
370+ const q = r + j + 1, hn = ch[rh + j + 1] * (du[q + 1] - du[q]);
371+ const fo = hn - hp + gvB[j] - gvA[j] - Kb * (dL1[q - P] + dL1[q + P] + dL1[q - 1] + dL1[q + 1] - 4 * dL1[q]);
372+ out[q] = fo; vv[q] += hh * fo; hp = hn;
373+ }
374+ gvA.set(gvB);
375+ }
376+ }
377+ // damp vin into vout, half kick with ff, drift uu
378+ function dampKickDrift(vin, vout, ff, uu) {
379+ for (let i = 1; i <= N; i++) {
380+ const r = i * P;
381+ for (let q = r + 1; q <= r + N; q++) {
382+ const x = vin[q], y = x - h * (c0 * x - c1 * (vin[q - P] + vin[q + P] + vin[q - 1] + vin[q + 1] - 4 * x)) + hh * ff[q];
383+ vout[q] = y; uu[q] += h * y;
384+ }
385+ }
386+ }
387+ function kickDrift(vv, ff, uu) {
388+ for (let i = 1; i <= N; i++) { const r = i * P; for (let q = r + 1; q <= r + N; q++) { const y = vv[q] + hh * ff[q]; vv[q] = y; uu[q] += h * y; } }
389+ }
390+ function damp(vin, vout) {
391+ for (let i = 1; i <= N; i++) {
392+ const r = i * P;
393+ for (let q = r + 1; q <= r + N; q++) { const x = vin[q]; vout[q] = x - h * (c0 * x - c1 * (vin[q - P] + vin[q + P] + vin[q - 1] + vin[q + 1] - 4 * x)); }
394+ }
395+ }
396+ function normalize(a, b) {
397+ let s = 0;
398+ for (let i = 1; i <= N; i++) for (let q = i * P + 1; q <= i * P + N; q++) s += a[q] * a[q] + b[q] * b[q];
399+ const n = Math.sqrt(s);
400+ for (let i = 1; i <= N; i++) for (let q = i * P + 1; q <= i * P + N; q++) { a[q] /= n; b[q] /= n; }
401+ return n;
402+ }
403+
404+ const nsub = Math.max(1, Math.ceil(p.dtOut / p.hmax)), h = p.dtOut / nsub, hh = 0.5 * h;
405+ const nOut = p.nOut, K = p.decim, nD = Math.floor(nOut / K), NB = p.bandOf.nb;
406+ const sig = new Float32Array(nOut);
407+ const snapN = opts.snap ? Math.min(nOut, p.snapN) : 0, snap = opts.snap ? new Float32Array(snapN * NN) : null;
408+ const bands = opts.modes ? new Float32Array(nD * NB) : null;
409+ const hi1 = new Float32Array(nD), hi4 = new Float32Array(nD), energy = new Float32Array(nD);
410+ const logsum = tan ? new Float64Array(nD) : null;
411+ const finalBands = new Float32Array(NB), es = new Float64Array(NN), bs = new Float64Array(NB), tmp = new Float64Array(NN);
412+ const Q = new Float64Array(NN), Pm = new Float64Array(NN);
413+ const pk = (p.lpi + 1) * P + p.lpj + 1;
414+ let L = 0;
415+ function modeEnergies() {
416+ // q = S w S, p = S v S (S symmetric)
417+ const S = sh.S;
418+ const proj = (x, out) => {
419+ for (let i = 0; i < N; i++) for (let b = 0; b < N; b++) { let s = 0; const r = (i + 1) * P + 1; for (let j = 0; j < N; j++) s += x[r + j] * S[j * N + b]; tmp[i * N + b] = s; }
420+ for (let a = 0; a < N; a++) for (let b = 0; b < N; b++) { let s = 0; for (let i = 0; i < N; i++) s += S[a * N + i] * tmp[i * N + b]; out[a * N + b] = s; }
421+ };
422+ proj(w, Q); proj(v, Pm);
423+ let tot = 0;
424+ for (let k = 0; k < NN; k++) { const e = 0.5 * (Pm[k] * Pm[k] + sh.om[k] * sh.om[k] * Q[k] * Q[k]); es[k] = e; tot += e; }
425+ bs.fill(0);
426+ for (let k = 0; k < NN; k++) bs[p.bandOf.idx[k]] += es[k] / (tot || 1);
427+ return tot;
428+ }
429+ const scratch = new Float64Array(PP);
430+ if (tan) { forceKickKeep(w, f, scratch, 0); tforceKick(dw, g, scratch, 0); } else forceKick(w, f, scratch, 0);
431+ for (let i = 0; i < nOut; i++) {
432+ for (let s = 0; s < nsub; s++) {
433+ // damping of the previous substep (that of the last substep of each output sample is applied below), kick, drift
434+ if (s) {
435+ dampKickDrift(v, vb, f, w); let t = v; v = vb; vb = t;
436+ if (tan) { dampKickDrift(dv, dvb, g, dw); t = dv; dv = dvb; dvb = t; }
437+ } else {
438+ kickDrift(v, f, w);
439+ if (tan) kickDrift(dv, g, dw);
440+ }
441+ if (tan) { forceKickKeep(w, f, v, hh); tforceKick(dw, g, dv, hh); } else forceKick(w, f, v, hh);
442+ }
443+ damp(v, vb); let t = v; v = vb; vb = t;
444+ if (tan) { damp(dv, dvb); t = dv; dv = dvb; dvb = t; }
445+ if (!Number.isFinite(v[P + 1]) || !Number.isFinite(v[pk])) return { blowup: true };
446+ if (tan) { L += Math.log(normalize(dw, dv)); tforceKick(dw, g, scratch, 0); }
447+ if (p.listenSum) { let sv = 0; for (let a = 1; a <= N; a++) for (let q = a * P + 1; q <= a * P + N; q++) sv += v[q]; sig[i] = sv; } else sig[i] = v[pk];
448+ if (i < snapN) for (let a = 0; a < N; a++) for (let b = 0; b < N; b++) snap[i * NN + a * N + b] = w[(a + 1) * P + b + 1];
449+ if ((i + 1) % K === 0) {
450+ const jd = (i + 1) / K - 1;
451+ if (jd < nD) {
452+ energy[jd] = modeEnergies();
453+ let a1 = 0, a4 = 0;
454+ for (let k = 0; k < NN; k++) { const fk = p.fk[k]; const sh_ = es[k] / (energy[jd] || 1); if (fk > 1000) a1 += sh_; if (fk > 4000) a4 += sh_; }
455+ hi1[jd] = a1; hi4[jd] = a4;
456+ if (bands) for (let b = 0; b < NB; b++) bands[jd * NB + b] = bs[b];
457+ if (logsum) logsum[jd] = L;
458+ }
459+ }
460+ }
461+ modeEnergies();
462+ for (let b = 0; b < NB; b++) finalBands[b] = bs[b];
463+ return { sig, snap, bands, hi1, hi4, energy, logsum, finalBands };
464+}
465+
466+onmessage = (ev) => {
467+ const { gen, p, k, lyap } = ev.data;
468+ const sh = setupSheet(p);
469+ if (lyap) {
470+ // strike A again, only for the Lyapunov exponent (run in parallel with the other strikes)
471+ let r = strikeRun(p, sh, 1, { tangent: true });
472+ if (r.blowup) r = strikeRun(Object.assign({}, p, { hmax: p.hmax / 3 }), sh, 1, { tangent: true });
473+ if (r.blowup) { postMessage({ gen, lyap: true, failed: true }); return; }
474+ postMessage({ gen, lyap: true, logsum: r.logsum }, [r.logsum.buffer]);
475+ return;
476+ }
477+ const opts = { snap: k === 0, modes: k === 0, tangent: false };
478+ let r = strikeRun(p, sh, 1 + k * p.delta, opts);
479+ if (r.blowup) r = strikeRun(Object.assign({}, p, { hmax: p.hmax / 3 }), sh, 1 + k * p.delta, opts);
480+ if (r.blowup) { postMessage({ gen, k, failed: true }); return; }
481+ const tr = [r.sig.buffer, r.hi1.buffer, r.hi4.buffer, r.energy.buffer, r.finalBands.buffer];
482+ r.logsum = null;
483+ if (r.snap) tr.push(r.snap.buffer);
484+ if (r.bands) tr.push(r.bands.buffer);
485+ postMessage({ gen, k, r }, tr);
486+};
487+</script>
488+
489+<script>
490+(() => {
491+ const SR = 44100, NS = 8;
492+ const $ = (id) => document.getElementById(id);
493+ const css = (name) => getComputedStyle(document.documentElement).getPropertyValue(name).trim();
494+ const fmt = (v, d = 2) => Number(v).toFixed(d);
495+
496+ const PRESETS = {
497+ soft: { E0: 1000, Kb: 25, Hd: 10, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 4, thi: 0.8, f1: 100, dur: 2.5 },
498+ gong: { E0: 10000, Kb: 25, Hd: 10, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 4, thi: 0.8, f1: 100, dur: 2.5 },
499+ crash: { E0: 30000, Kb: 25, Hd: 10, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 4, thi: 0.8, f1: 100, dur: 2.5 },
500+ flat: { E0: 10000, Kb: 25, Hd: 0, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 4, thi: 0.8, f1: 100, dur: 2.5 },
501+ // chaotic: strikes differing by 1e-6 end with different spectra (scan of 2026-10-05, chaos_scan.py)
502+ cgong: { E0: 30000, Kb: 8, Hd: 10, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 6, thi: 2, f1: 150, dur: 2.5 },
503+ ccrash: { E0: 100000, Kb: 25, Hd: 10, N: 24, sx: 0.35, sy: 0.42, rad: 1.6, dl: -6, tlo: 6, thi: 2, f1: 100, dur: 2.5 },
504+ };
505+ const eFromSlider = (v) => Math.pow(10, 2 + 3 * v / 1000);
506+ const sliderFromE = (E) => Math.round((Math.log10(E) - 2) / 3 * 1000);
507+
508+ function sheetFreqs(N, Kb) {
509+ const M = N + 1, s = [];
510+ for (let k = 0; k < N; k++) s.push(4 * Math.sin((k + 1) * Math.PI / (2 * M)) ** 2);
511+ const om = new Float64Array(N * N);
512+ for (let a = 0; a < N; a++) for (let b = 0; b < N; b++) { const lam = s[a] + s[b]; om[a * N + b] = Math.sqrt(lam + Kb * lam * lam); }
513+ return om;
514+ }
515+
516+ function readParams() {
517+ const E0 = eFromSlider(+$("E0").value), Kb = +$("Kb").value, Hd = +$("Hd").value, N = +$("N").value;
518+ const sx = +$("sx").value, sy = +$("sy").value, rad = +$("rad").value, dl = +$("dl").value;
519+ const tlo = +$("tlo").value, thi = Math.min(+$("thi").value, tlo), dur = +$("dur").value;
520+ const om = sheetFreqs(N, Kb);
521+ let om1 = Infinity, omN = 0;
522+ for (const o of om) { om1 = Math.min(om1, o); omN = Math.max(omN, o); }
523+ const f1max = 19000 / (omN / om1), f1 = Math.min(+$("f1").value, f1max), ftop = f1 * omN / om1;
524+ const M = N + 1, ci = Math.max(1, Math.min(N, Math.round(sx * M))), cj = Math.max(1, Math.min(N, Math.round(sy * M)));
525+ const lpi = Math.max(0, Math.min(N - 1, Math.round(+$("lx").value * M) - 1)), lpj = Math.max(0, Math.min(N - 1, Math.round(+$("ly").value * M) - 1));
526+ const sum = $("lisSum").getAttribute("aria-pressed") === "true";
527+ $("E0v").textContent = E0 >= 1000 ? fmt(E0 / 1000, E0 >= 10000 ? 0 : 1) + "k" : fmt(E0, 0);
528+ $("Kbv").textContent = fmt(Kb, 1);
529+ $("Hdv").textContent = fmt(Hd, 1);
530+ $("Nv").textContent = N + " × " + N;
531+ $("sposv").textContent = "bead (" + ci + ", " + cj + ")";
532+ $("radv").textContent = fmt(rad, 1) + " beads";
533+ $("dlv").textContent = "1e" + dl + " in strength";
534+ $("tlov").textContent = fmt(tlo, 1) + " s";
535+ $("thiv").textContent = fmt(thi, 2) + " s";
536+ $("f1v").textContent = Math.round(f1) + " Hz · top " + (ftop >= 1000 ? fmt(ftop / 1000, 1) + "k" : Math.round(ftop));
537+ $("f1note").hidden = !(+$("f1").value > f1max);
538+ $("durv").textContent = fmt(dur, 1) + " s";
539+ $("lisv").textContent = sum ? "whole sheet" : "bead (" + (lpi + 1) + ", " + (lpj + 1) + ")";
540+ $("lx").disabled = sum; $("ly").disabled = sum;
541+ const tau = om1 / (2 * Math.PI * f1);
542+ // damping: mode decay (amplitude) about (c0 + c1 lam) / 2, lam the Laplacian eigenvalue, from the ring times
543+ const s1 = 4 * Math.sin(Math.PI / (2 * M)) ** 2, sN = 4 * Math.sin(N * Math.PI / (2 * M)) ** 2, lam1 = 2 * s1, lamN = 2 * sN;
544+ const sLo = Math.log(1000) / tlo * tau, sHi = Math.log(1000) / thi * tau;
545+ const c1 = Math.max(0, 2 * (sHi - sLo) / (lamN - lam1)), c0 = Math.max(0, 2 * sLo - c1 * lam1);
546+ const fk = Array.from(om, (o) => f1 * o / om1);
547+ // log-spaced bands from f1 to the top
548+ const NB = 36, idx = new Int32Array(N * N), edges = [];
549+ for (let b = 0; b <= NB; b++) edges.push(f1 * Math.pow(ftop * 1.0001 / f1, b / NB));
550+ for (let k = 0; k < N * N; k++) { let b = Math.floor(Math.log(fk[k] / f1) / Math.log(ftop * 1.0001 / f1) * NB); idx[k] = Math.max(0, Math.min(NB - 1, b)); }
551+ return {
552+ E0, Kb, Hd, N, T: 1, beta: 1, ci, cj, rad, delta: Math.pow(10, dl), tau, c0, c1, dur, f1, ftop, fk,
553+ bandOf: { nb: NB, idx, edges }, lpi, lpj, listenSum: sum,
554+ dtOut: 1 / SR / tau, nOut: Math.round(dur * SR), decim: Math.round(SR / 400),
555+ // leapfrog at about 0.6 of its linear stability limit; the band statistics agree with 4x smaller steps
556+ // (checked 2026-10-05), and a strike that blows up is rerun with a third of the step
557+ hmax: 1.2 / omN, snapN: Math.round(0.08 * SR),
558+ };
559+ }
560+ function setControls(p) {
561+ $("E0").value = sliderFromE(p.E0); $("Kb").value = p.Kb; $("Hd").value = p.Hd; $("N").value = p.N; $("sx").value = p.sx; $("sy").value = p.sy;
562+ $("rad").value = p.rad; $("dl").value = p.dl; $("tlo").value = p.tlo; $("thi").value = p.thi; $("f1").value = p.f1; $("dur").value = p.dur;
563+ }
564+
565+ // ---------------- workers ----------------
566+ const workerURL = URL.createObjectURL(new Blob([$("worker-src").textContent], { type: "text/javascript" }));
567+ const NW = Math.max(2, Math.min(NS + 1, (navigator.hardwareConcurrency || 4) - 1));
568+ let workers = [], gen = 0, params = null, S = [], queue = [], done = 0, t0run = 0, lamA = null;
569+ function run() {
570+ stopAudio(); stopSlow(true);
571+ params = readParams();
572+ gen++;
573+ for (const w of workers) w.terminate();
574+ workers = []; S = new Array(NS).fill(null); done = 0; bufs = []; lamA = null;
575+ queue = ["lyap"].concat(Array.from({ length: NS }, (_, k) => k));
576+ t0run = performance.now();
577+ for (let i = 0; i < Math.min(NW, NS); i++) { const w = new Worker(workerURL); w.onmessage = onResult; workers.push(w); next(w); }
578+ $("status").textContent = "Computing 8 strikes in parallel…";
579+ renderSelector(); updateStats(); drawAll();
580+ }
581+ function next(w) {
582+ if (!queue.length) return;
583+ const k = queue.shift();
584+ if (k === "lyap") w.postMessage({ gen, p: params, k: 0, lyap: true }); else w.postMessage({ gen, p: params, k });
585+ }
586+ function onResult(ev) {
587+ const m = ev.data;
588+ if (m.gen !== gen) return;
589+ next(ev.target);
590+ if (m.failed) { $("status").textContent = "The strike blew up numerically; try a softer strike."; return; }
591+ if (m.lyap) {
592+ lamA = localLyap(m.logsum);
593+ if (S[0]) S[0].lam = lamA;
594+ done++;
595+ } else {
596+ const r = m.r; S[m.k] = r; done++;
597+ if (m.k === 0) { r.lam = lamA; audioGain = 0.8 / (maxAbs(r.sig) || 1); }
598+ }
599+ if (S[0]) for (let k = 0; k < NS; k++) if (S[k] && !S[k].audio) S[k].audio = scaleAudio(S[k].sig);
600+ $("status").textContent = done < NS + 1 ? "Computed " + done + " of 9 runs (8 strikes and the separation rate)…" : "Ready (" + fmt((performance.now() - t0run) / 1000, 1) + " s). Strikes differ by 1e" + $("dl").value + " in strength.";
601+ renderSelector(); updateStats(); drawAll();
602+ }
603+ const maxAbs = (a) => { let m = 0; for (let i = 0; i < a.length; i++) m = Math.max(m, Math.abs(a[i])); return m; };
604+ function localLyap(logsum) {
605+ const dtD = params.decim * params.dtOut, w = Math.max(1, Math.round(0.05 / (params.decim / SR)));
606+ const out = new Float32Array(logsum.length);
607+ for (let i = 0; i < logsum.length; i++) { const j = Math.max(0, i - w); out[i] = i > j ? (logsum[i] - logsum[j]) / ((i - j) * dtD) / params.tau : NaN; }
608+ return out;
609+ }
610+ function chaosEnd(lam) { let last = -1; for (let i = 0; i < lam.length; i++) if (lam[i] > 40) last = i; return last < 0 ? null : (last + 1) * params.decim / SR; }
611+ function updateStats() {
612+ const A = S[0];
613+ if (!A) { for (const id of ["sBloom", "sBloomT", "sChaos", "sLam", "sSpread"]) $(id).textContent = "–"; return; }
614+ let pk = 0, pi = 0;
615+ for (let i = 0; i < A.hi4.length; i++) if (A.hi4[i] > pk) { pk = A.hi4[i]; pi = i; }
616+ $("sBloom").textContent = fmt(100 * pk, pk < 0.01 ? 2 : 0) + "%";
617+ $("sBloomT").textContent = fmt((pi + 1) * params.decim / SR * 1000, 0) + " ms";
618+ if (!A.lam) { $("sChaos").textContent = "computing…"; $("sLam").textContent = "computing…"; }
619+ const end = A.lam ? chaosEnd(A.lam) : null;
620+ if (A.lam) $("sChaos").textContent = end === null ? "never chaotic" : end > params.dur - 0.15 ? "still chaotic at the end" : fmt(end, 2) + " s";
621+ let peak = 0; if (A.lam) for (const v of A.lam) if (v > peak) peak = v;
622+ if (A.lam) $("sLam").textContent = peak > 40 ? fmt(peak, 0) + " /s" : "≈ 0";
623+ if (S.some((x) => !x)) { $("sSpread").textContent = "computing…"; return; }
624+ const NB = params.bandOf.nb, sd = [];
625+ for (let b = 0; b < NB; b++) {
626+ const mean = S.reduce((s, x) => s + x.finalBands[b], 0) / NS;
627+ if (mean < 1e-3) continue;
628+ const vals = S.map((x) => Math.log10(Math.max(x.finalBands[b], 1e-12))), m = vals.reduce((s, v) => s + v, 0) / NS;
629+ sd.push(Math.sqrt(vals.reduce((s, v) => s + (v - m) ** 2, 0) / NS));
630+ }
631+ sd.sort((a, b) => a - b);
632+ const med = sd.length ? sd[Math.floor(sd.length / 2)] : 0;
633+ $("sSpread").textContent = med < 0.01 ? "identical" : "differ by " + fmt(10 * med, 1) + " dB";
634+ }
635+
636+ // ---------------- selector ----------------
637+ let sel = 0;
638+ const strikeColor = (k) => (k === 0 ? css("--x") : k === 1 ? css("--y") : css("--ink"));
639+ const strikeName = (k) => (k === 0 ? "A" : k === 1 ? "B" : String(k + 1));
640+ function renderSelector() {
641+ const box = $("strikeSel");
642+ if (!box.childElementCount) {
643+ for (let k = 0; k < NS; k++) { const b = document.createElement("button"); b.type = "button"; b.className = "sk"; b.id = "sk" + k; b.innerHTML = '<span class="sw"></span>' + strikeName(k); b.addEventListener("click", () => selectStrike(k)); box.appendChild(b); }
644+ const pb = document.createElement("button"); pb.type = "button"; pb.className = "btn"; pb.id = "playSel"; pb.textContent = "Play this strike"; pb.addEventListener("click", () => play(sel, null)); box.appendChild(pb);
645+ }
646+ for (let k = 0; k < NS; k++) { const b = $("sk" + k); b.disabled = !S[k]; b.setAttribute("aria-pressed", String(k === sel)); b.querySelector(".sw").style.background = strikeColor(k); }
647+ $("playSel").disabled = !S[sel];
648+ }
649+ function selectStrike(k) { if (!S[k]) return; sel = k; renderSelector(); drawHi(); drawFinal(); drawSpec(); }
650+ document.addEventListener("keydown", (e) => {
651+ if (e.target && (e.target.tagName === "INPUT" || e.target.tagName === "TEXTAREA")) return;
652+ if (e.key !== "ArrowRight" && e.key !== "ArrowLeft") return;
653+ const dir = e.key === "ArrowRight" ? 1 : -1;
654+ for (let i = 1; i <= NS; i++) { const k = (sel + dir * i + NS) % NS; if (S[k]) { selectStrike(k); e.preventDefault(); return; } }
655+ });
656+
657+ // ---------------- canvas helpers ----------------
658+ function setup(canvas) {
659+ const r = canvas.getBoundingClientRect(), dpr = window.devicePixelRatio || 1;
660+ const w = Math.max(1, Math.round(r.width * dpr)), h = Math.max(1, Math.round(r.height * dpr));
661+ if (canvas.width !== w || canvas.height !== h) { canvas.width = w; canvas.height = h; }
662+ const ctx = canvas.getContext("2d"); ctx.setTransform(dpr, 0, 0, dpr, 0, 0); ctx.clearRect(0, 0, r.width, r.height);
663+ return { ctx, W: r.width, H: r.height };
664+ }
665+ const PAD = { l: 48, r: 12, t: 10, b: 26 };
666+ function niceStep(v) { const p = Math.pow(10, Math.floor(Math.log10(v))), m = v / p; return (m < 1.5 ? 1 : m < 3.5 ? 2 : m < 7.5 ? 5 : 10) * p; }
667+ const hz = (f) => (f >= 1000 ? fmt(f / 1000, f >= 10000 || f % 1000 === 0 ? 0 : 1) + "k" : String(Math.round(f)));
668+ function axesTime(ctx, W, H, ymin, ymax, yticks, yfmt, dur) {
669+ const x0 = PAD.l, x1 = W - PAD.r, y0 = PAD.t, y1 = H - PAD.b;
670+ ctx.strokeStyle = css("--grid"); ctx.lineWidth = 1; ctx.fillStyle = css("--muted"); ctx.font = "11px " + css("--mono");
671+ ctx.textAlign = "right"; ctx.textBaseline = "middle";
672+ for (const v of yticks) { const y = y1 - (v - ymin) / (ymax - ymin) * (y1 - y0); ctx.beginPath(); ctx.moveTo(x0, y); ctx.lineTo(x1, y); ctx.stroke(); ctx.fillText(yfmt(v), x0 - 6, y); }
673+ ctx.textAlign = "center"; ctx.textBaseline = "top";
674+ const step = dur <= 2 ? 0.25 : 0.5;
675+ for (let t = 0; t <= dur + 1e-9; t += step) { const x = x0 + t / dur * (x1 - x0); ctx.beginPath(); ctx.moveTo(x, y0); ctx.lineTo(x, y1); ctx.stroke(); ctx.fillText(fmt(t, step < 0.5 ? 2 : 1) + (t === 0 ? " s" : ""), x, y1 + 6); }
676+ return { x0, x1, y0, y1, X: (t) => x0 + t / dur * (x1 - x0), Y: (v) => y1 - (v - ymin) / (ymax - ymin) * (y1 - y0) };
677+ }
678+ function chaosSpan(ctx, ax) { if (!S[0] || !S[0].lam) return; const end = chaosEnd(S[0].lam); if (end === null) return; ctx.fillStyle = css("--chaos-wash"); ctx.fillRect(ax.x0, ax.y0, ax.X(end) - ax.x0, ax.y1 - ax.y0); }
679+ function seriesTime(ctx, ax, arr, color, width, dash) {
680+ if (!arr) return;
681+ const dt = params.decim / SR;
682+ ctx.strokeStyle = color; ctx.lineWidth = width; ctx.setLineDash(dash || []); ctx.beginPath();
683+ let started = false;
684+ for (let i = 0; i < arr.length; i++) { const v = arr[i]; if (!Number.isFinite(v)) { started = false; continue; } const x = ax.X((i + 1) * dt), y = ax.Y(v); if (!started) { ctx.moveTo(x, y); started = true; } else ctx.lineTo(x, y); }
685+ ctx.stroke(); ctx.setLineDash([]);
686+ }
687+ function hexToRgb(h) { h = h.replace("#", ""); if (h.length === 3) h = h.split("").map((c) => c + c).join(""); const n = parseInt(h, 16); return [(n >> 16) & 255, (n >> 8) & 255, n & 255]; }
688+ function mix(a, b, t) { return [a[0] + (b[0] - a[0]) * t, a[1] + (b[1] - a[1]) * t, a[2] + (b[2] - a[2]) * t]; }
689+ function heatColor(v, bg, ink, hot) { return v < 0.72 ? mix(bg, ink, v / 0.72) : mix(ink, hot, (v - 0.72) / 0.28); }
690+
691+ // ---------------- the sheet ----------------
692+ let smOn = false, smIdx = 0, smRaf = 0, smLast = 0;
693+ function sheetGeom(W, H) { const N = params.N, pad = 14, side = Math.min(W, H) - 2 * pad, cell = side / (N + 1); return { N, pad, side, cell, x0: (W - side) / 2, y0: (H - side) / 2 }; }
694+ function drawSheet() {
695+ const { ctx, W, H } = setup($("cSheet"));
696+ if (!params) return;
697+ const g = sheetGeom(W, H), N = g.N, A = S[0];
698+ // frame (pinned edge)
699+ ctx.fillStyle = css("--rule"); ctx.fillRect(g.x0, g.y0, g.side, g.side);
700+ const bg = hexToRgb(css("--panel")), up = hexToRgb(css("--up")), dn = hexToRgb(css("--down"));
701+ let wm = 0;
702+ const w = new Float32Array(N * N);
703+ if (A && A.snap) {
704+ const nS = A.snap.length / (N * N);
705+ for (let i = 0; i < A.snap.length; i += 7) wm = Math.max(wm, Math.abs(A.snap[i]));
706+ if (smOn) { const i = Math.min(nS - 1, Math.floor(smIdx)); for (let k = 0; k < N * N; k++) w[k] = A.snap[i * N * N + k]; }
707+ }
708+ // rest dome: faint concentric shading
709+ for (let i = 0; i < N; i++) for (let j = 0; j < N; j++) {
710+ const X = g.x0 + (j + 0.5 + 0.5) * g.cell, Y = g.y0 + (i + 0.5 + 0.5) * g.cell;
711+ const t = wm > 0 ? w[i * N + j] / wm : 0;
712+ const c = t >= 0 ? mix(bg, up, Math.min(1, t * 1.6)) : mix(bg, dn, Math.min(1, -t * 1.6));
713+ ctx.fillStyle = "rgb(" + c.map(Math.round).join(",") + ")";
714+ ctx.fillRect(X - g.cell * 0.46, Y - g.cell * 0.46, g.cell * 0.92, g.cell * 0.92);
715+ }
716+ if (params.Hd > 0) {
717+ ctx.strokeStyle = css("--muted"); ctx.globalAlpha = 0.35; ctx.lineWidth = 1;
718+ for (const fr of [0.25, 0.5, 0.75]) {
719+ // contour of the dome x(1-x)y(1-y) = fr / 16: draw by marching a polar parametrization
720+ ctx.beginPath();
721+ for (let a = 0; a <= 64; a++) {
722+ const th = 2 * Math.PI * a / 64; let lo = 0, hi = 0.5;
723+ for (let it = 0; it < 30; it++) { const r = (lo + hi) / 2, x = 0.5 + r * Math.cos(th), y = 0.5 + r * Math.sin(th); if (16 * x * (1 - x) * y * (1 - y) > fr) lo = r; else hi = r; }
724+ const r = (lo + hi) / 2, X = g.x0 + g.cell * 0.5 + (0.5 + r * Math.cos(th)) * g.cell * (N + 1) - g.cell * 0.5, Y = g.y0 + (0.5 + r * Math.sin(th)) * g.cell * (N + 1);
725+ if (a === 0) ctx.moveTo(X, Y); else ctx.lineTo(X, Y);
726+ }
727+ ctx.stroke();
728+ }
729+ ctx.globalAlpha = 1;
730+ }
731+ // strike patch
732+ const sX = g.x0 + (params.cj + 0.5) * g.cell, sY = g.y0 + (params.ci + 0.5) * g.cell;
733+ if (!smOn) { ctx.strokeStyle = css("--chaos"); ctx.lineWidth = 3; ctx.beginPath(); ctx.arc(sX, sY, Math.max(6, params.rad * g.cell * 1.2), 0, 2 * Math.PI); ctx.stroke(); }
734+ // listener
735+ ctx.strokeStyle = css("--y"); ctx.lineWidth = 2.2;
736+ if (params.listenSum) { ctx.strokeRect(g.x0 + 3, g.y0 + 3, g.side - 6, g.side - 6); }
737+ else { const lX = g.x0 + (params.lpj + 1.5) * g.cell, lY = g.y0 + (params.lpi + 1.5) * g.cell; ctx.beginPath(); ctx.arc(lX, lY, Math.max(6, g.cell * 0.8), 0, 2 * Math.PI); ctx.stroke(); }
738+ ctx.fillStyle = css("--ink"); ctx.font = "11px " + css("--mono"); ctx.textAlign = "left"; ctx.textBaseline = "top";
739+ if (smOn && A) ctx.fillText("t = " + fmt(smIdx / SR * 1000, 2) + " ms", g.x0 + 6, g.y0 + 4);
740+ }
741+ function smTick(now) {
742+ const A = S[0];
743+ if (!smOn || !A) return;
744+ const dtReal = smLast ? Math.min(0.1, (now - smLast) / 1000) : 0;
745+ smLast = now; smIdx += 0.004 * dtReal * SR;
746+ if (smIdx >= A.snap.length / (params.N * params.N) - 1) { stopSlow(); return; }
747+ drawSheet(); smRaf = requestAnimationFrame(smTick);
748+ }
749+ function startSlow() { if (!S[0]) return; smOn = true; smIdx = 0; smLast = 0; $("slow").setAttribute("aria-pressed", "true"); $("slow").textContent = "Stop slow motion"; smRaf = requestAnimationFrame(smTick); }
750+ function stopSlow(silent) { smOn = false; cancelAnimationFrame(smRaf); $("slow").setAttribute("aria-pressed", "false"); $("slow").textContent = "Slow motion"; if (!silent) drawSheet(); }
751+ $("cSheet").addEventListener("click", (e) => {
752+ if (!params) return;
753+ const r = $("cSheet").getBoundingClientRect(), g = sheetGeom(r.width, r.height);
754+ const fx = (e.clientX - r.left - g.x0) / g.side, fy = (e.clientY - r.top - g.y0) / g.side;
755+ if (fx < 0 || fx > 1 || fy < 0 || fy > 1) return;
756+ if (e.shiftKey) {
757+ $("lx").value = Math.min(0.95, Math.max(0.05, fy)); $("ly").value = Math.min(0.95, Math.max(0.05, fx));
758+ if ($("lisSum").getAttribute("aria-pressed") === "true") { $("lisSum").setAttribute("aria-pressed", "false"); $("lisOne").setAttribute("aria-pressed", "true"); }
759+ } else { $("sx").value = Math.min(0.92, Math.max(0.08, fy)); $("sy").value = Math.min(0.92, Math.max(0.08, fx)); }
760+ readParams(); run();
761+ });
762+
763+ // ---------------- plots ----------------
764+ function drawModes() {
765+ const { ctx, W, H } = setup($("cModes"));
766+ const ax = axesTime(ctx, W, H, 0, 1, [], () => "", params.dur);
767+ const A = S[0];
768+ if (!A || !A.bands) return;
769+ const NB = params.bandOf.nb, nD = A.hi4.length, off = document.createElement("canvas");
770+ off.width = nD; off.height = NB;
771+ const octx = off.getContext("2d"), img = octx.createImageData(nD, NB);
772+ const bg = hexToRgb(css("--panel")), ink = hexToRgb(css("--ink")), hot = hexToRgb(css("--chaos"));
773+ for (let i = 0; i < nD; i++) for (let b = 0; b < NB; b++) {
774+ const v = Math.min(1, Math.max(0, (Math.log10(Math.max(A.bands[i * NB + b], 1e-12)) + 6) / 6));
775+ const c = heatColor(v, bg, ink, hot), o = ((NB - 1 - b) * nD + i) * 4;
776+ img.data[o] = c[0]; img.data[o + 1] = c[1]; img.data[o + 2] = c[2]; img.data[o + 3] = 255;
777+ }
778+ octx.putImageData(img, 0, 0);
779+ ctx.imageSmoothingEnabled = false;
780+ ctx.drawImage(off, ax.x0, ax.y0, ax.x1 - ax.x0, ax.y1 - ax.y0);
781+ ctx.fillStyle = css("--muted"); ctx.font = "11px " + css("--mono"); ctx.textAlign = "right"; ctx.textBaseline = "middle";
782+ const lf0 = Math.log(params.f1), lf1 = Math.log(params.ftop);
783+ for (const f of [100, 200, 500, 1000, 2000, 4000, 8000, 16000]) if (f >= params.f1 && f <= params.ftop) ctx.fillText(hz(f), ax.x0 - 6, ax.y1 - (Math.log(f) - lf0) / (lf1 - lf0) * (ax.y1 - ax.y0));
784+ }
785+ function drawHi() {
786+ const { ctx, W, H } = setup($("cHi"));
787+ const ax = axesTime(ctx, W, H, 0, 1, [0, 0.25, 0.5, 0.75, 1], (v) => fmt(v * 100, 0) + "%", params.dur);
788+ chaosSpan(ctx, ax);
789+ for (let k = NS - 1; k >= 0; k--) { const x = S[k]; if (!x) continue; const col = k === 0 ? css("--x") : k === 1 ? css("--y") : css("--faint"); seriesTime(ctx, ax, x.hi1, col, 0.8); seriesTime(ctx, ax, x.hi4, col, k <= 1 ? 1.6 : 1); }
790+ if (S[sel]) { seriesTime(ctx, ax, S[sel].hi4, strikeColor(sel), 3); seriesTime(ctx, ax, S[sel].hi1, strikeColor(sel), 1.4, [4, 3]); }
791+ }
792+ function drawFinal() {
793+ const { ctx, W, H } = setup($("cFinal"));
794+ if (!params) return;
795+ const NB = params.bandOf.nb, e = params.bandOf.edges, x0 = PAD.l, x1 = W - PAD.r, y0 = PAD.t, y1 = H - PAD.b, lo = -8, hi = 0;
796+ const lf0 = Math.log10(e[0]), lf1 = Math.log10(e[NB]);
797+ const X = (f) => x0 + (Math.log10(f) - lf0) / (lf1 - lf0) * (x1 - x0), Y = (v) => y1 - (v - lo) / (hi - lo) * (y1 - y0);
798+ ctx.strokeStyle = css("--grid"); ctx.lineWidth = 1; ctx.fillStyle = css("--muted"); ctx.font = "11px " + css("--mono");
799+ ctx.textAlign = "right"; ctx.textBaseline = "middle";
800+ for (let v = lo; v <= hi; v += 2) { ctx.beginPath(); ctx.moveTo(x0, Y(v)); ctx.lineTo(x1, Y(v)); ctx.stroke(); ctx.fillText("1e" + v, x0 - 6, Y(v)); }
801+ ctx.textAlign = "center"; ctx.textBaseline = "top";
802+ for (const f of [100, 200, 500, 1000, 2000, 5000, 10000]) if (f >= e[0] && f <= e[NB]) { ctx.beginPath(); ctx.moveTo(X(f), y0); ctx.lineTo(X(f), y1); ctx.stroke(); ctx.fillText(hz(f), X(f), y1 + 6); }
803+ const line = (s, color, w) => { if (!s) return; ctx.strokeStyle = color; ctx.lineWidth = w; ctx.beginPath(); for (let b = 0; b < NB; b++) { const fc = Math.sqrt(e[b] * e[b + 1]), v = Math.max(lo, Math.log10(Math.max(s[b], 1e-12))); if (b === 0) ctx.moveTo(X(fc), Y(v)); else ctx.lineTo(X(fc), Y(v)); } ctx.stroke(); };
804+ for (let k = 2; k < NS; k++) if (S[k]) line(S[k].finalBands, css("--faint"), 1);
805+ if (S[1]) line(S[1].finalBands, css("--y"), 1.4);
806+ if (S[0]) line(S[0].finalBands, css("--x"), 1.6);
807+ if (S[sel]) line(S[sel].finalBands, strikeColor(sel), 3);
808+ }
809+ function drawLam() {
810+ const { ctx, W, H } = setup($("cLam"));
811+ const A = S[0];
812+ if (!A || !A.lam) { axesTime(ctx, W, H, 0, 1, [0, 1], () => "", params.dur); return; }
813+ let mx = 0; for (const v of A.lam) if (Number.isFinite(v) && v > mx) mx = v;
814+ mx = Math.max(mx * 1.1, 60);
815+ const st = niceStep(mx / 3), ticks = []; for (let v = 0; v <= mx; v += st) ticks.push(v);
816+ const ax = axesTime(ctx, W, H, -0.15 * mx, mx, ticks, (v) => fmt(v, 0), params.dur);
817+ chaosSpan(ctx, ax);
818+ ctx.strokeStyle = css("--rule"); ctx.beginPath(); ctx.moveTo(ax.x0, ax.Y(0)); ctx.lineTo(ax.x1, ax.Y(0)); ctx.stroke();
819+ const e = A.energy, e0 = e[0] || 1, eDb = new Float32Array(e.length);
820+ for (let i = 0; i < e.length; i++) eDb[i] = -0.15 * mx + (Math.max(10 * Math.log10(e[i] / e0), -80) + 80) / 80 * 1.15 * mx;
821+ seriesTime(ctx, ax, eDb, css("--muted"), 1, [4, 4]);
822+ seriesTime(ctx, ax, A.lam, css("--x"), 1.6);
823+ }
824+
825+ // ---------------- spectrogram ----------------
826+ function fft(re, im) {
827+ const n = re.length;
828+ for (let i = 1, j = 0; i < n; i++) { let bit = n >> 1; for (; j & bit; bit >>= 1) j ^= bit; j ^= bit; if (i < j) { [re[i], re[j]] = [re[j], re[i]]; [im[i], im[j]] = [im[j], im[i]]; } }
829+ for (let len = 2; len <= n; len <<= 1) { const ang = -2 * Math.PI / len, wr = Math.cos(ang), wi = Math.sin(ang); for (let i = 0; i < n; i += len) { let cr = 1, ci = 0; for (let k = 0; k < len / 2; k++) { const a = i + k, b = a + len / 2; const tr = re[b] * cr - im[b] * ci, ti = re[b] * ci + im[b] * cr; re[b] = re[a] - tr; im[b] = im[a] - ti; re[a] += tr; im[a] += ti; const ncr = cr * wr - ci * wi; ci = cr * wi + ci * wr; cr = ncr; } } }
830+ }
831+ const specFmax = () => Math.min(SR / 2, Math.max(2000, 1.1 * params.ftop));
832+ function spectrogramOf(x) {
833+ if (x.spec) return x.spec;
834+ const sig = x.sig, N = 1024, hop = 256, frames = Math.max(1, Math.floor((sig.length - N) / hop)), bins = Math.floor(specFmax() / (SR / N));
835+ const win = new Float32Array(N); for (let i = 0; i < N; i++) win[i] = 0.5 - 0.5 * Math.cos(2 * Math.PI * i / (N - 1));
836+ const mag = new Float32Array(frames * bins), re = new Float64Array(N), im = new Float64Array(N);
837+ let mx = -Infinity;
838+ for (let f = 0; f < frames; f++) { for (let i = 0; i < N; i++) { re[i] = sig[f * hop + i] * win[i]; im[i] = 0; } fft(re, im); for (let b = 0; b < bins; b++) { const v = 10 * Math.log10(re[b] * re[b] + im[b] * im[b] + 1e-20); mag[f * bins + b] = v; if (v > mx) mx = v; } }
839+ x.spec = { mag, frames, bins, mx };
840+ return x.spec;
841+ }
842+ function drawSpec() {
843+ const { ctx, W, H } = setup($("cSpec"));
844+ const ax = axesTime(ctx, W, H, 0, 1, [], () => "", params.dur);
845+ const X = S[sel] || S[0];
846+ if (!X || !S[0]) return;
847+ const ref = spectrogramOf(S[0]).mx, { mag, frames, bins } = spectrogramOf(X);
848+ const off = document.createElement("canvas"); off.width = frames; off.height = bins;
849+ const octx = off.getContext("2d"), img = octx.createImageData(frames, bins);
850+ const bg = hexToRgb(css("--panel")), ink = hexToRgb(css("--ink")), hot = hexToRgb(css("--chaos"));
851+ for (let f = 0; f < frames; f++) for (let b = 0; b < bins; b++) { const v = Math.min(1, Math.max(0, (mag[f * bins + b] - (ref - 80)) / 80)); const c = heatColor(v, bg, ink, hot), o = ((bins - 1 - b) * frames + f) * 4; img.data[o] = c[0]; img.data[o + 1] = c[1]; img.data[o + 2] = c[2]; img.data[o + 3] = 255; }
852+ octx.putImageData(img, 0, 0);
853+ const t0 = 512 / SR, t1 = (frames * 256 + 512) / SR;
854+ ctx.drawImage(off, ax.X(t0), ax.y0, ax.X(t1) - ax.X(t0), ax.y1 - ax.y0);
855+ ctx.fillStyle = css("--muted"); ctx.font = "11px " + css("--mono"); ctx.textAlign = "right"; ctx.textBaseline = "middle";
856+ const fm = specFmax(), step = niceStep(fm / 4);
857+ for (let f = 0; f <= fm; f += step) ctx.fillText(hz(f), ax.x0 - 6, ax.y1 - f / fm * (ax.y1 - ax.y0));
858+ const ks = S[sel] ? sel : 0;
859+ ctx.fillStyle = strikeColor(ks); ctx.font = "600 12px " + css("--sans"); ctx.textAlign = "right"; ctx.textBaseline = "top"; ctx.fillText("strike " + strikeName(ks), ax.x1 - 6, ax.y0 + 4);
860+ }
861+ function drawAll() { if (!params) return; drawSheet(); drawModes(); drawHi(); drawFinal(); drawLam(); drawSpec(); drawOverlays(null, null); }
862+
863+ // ---------------- audio ----------------
864+ let actx = null, bufs = [], src = null, playStart = 0, playing = null, raf = 0, audioGain = 1;
865+ function scaleAudio(s) { const o = new Float32Array(s.length), fade = Math.round(0.01 * SR); for (let i = 0; i < s.length; i++) { const tail = s.length - 1 - i; o[i] = Math.max(-1, Math.min(1, s[i] * audioGain * (tail < fade ? tail / fade : 1))); } return o; }
866+ function ensureCtx() { if (!actx) actx = new (window.AudioContext || window.webkitAudioContext)(); if (actx.state === "suspended") actx.resume(); }
867+ function toBuffer(d) { const b = actx.createBuffer(1, d.length, SR); b.copyToChannel(d, 0); return b; }
868+ function stopAudio() { if (src) { try { src.onended = null; src.stop(); } catch (e) {} src = null; } playing = null; cancelAnimationFrame(raf); for (const id of ["playA", "playB", "playAB", "playSel"]) if ($(id)) $(id).setAttribute("aria-pressed", "false"); drawOverlays(null, null); }
869+ function play(k, then) {
870+ const x = S[k]; if (!x || !x.audio) return;
871+ ensureCtx(); stopAudio();
872+ if (!bufs[k]) bufs[k] = toBuffer(x.audio);
873+ src = actx.createBufferSource(); src.buffer = bufs[k]; src.connect(actx.destination);
874+ const nx = then === null || then === undefined ? null : then;
875+ src.onended = () => { src = null; if (nx !== null) setTimeout(() => play(nx, null), 350); else stopAudio(); };
876+ playStart = actx.currentTime; playing = k; src.start();
877+ if (k === 0) $("playA").setAttribute("aria-pressed", "true"); else if (k === 1) $("playB").setAttribute("aria-pressed", "true");
878+ if (k === sel && $("playSel")) $("playSel").setAttribute("aria-pressed", "true");
879+ if (nx !== null) $("playAB").setAttribute("aria-pressed", "true");
880+ tick();
881+ }
882+ function tick() { if (playing === null) return; drawOverlays(actx.currentTime - playStart, playing); raf = requestAnimationFrame(tick); }
883+ function drawOverlays(t, which) {
884+ for (const [id, base] of [["cModesO", "cModes"], ["cHiO", "cHi"], ["cLamO", "cLam"], ["cSpecO", "cSpec"]]) {
885+ const { ctx, W } = setup($(id)); if (t === null) continue;
886+ const H = $(base).getBoundingClientRect().height, x = PAD.l + Math.min(t, params.dur) / params.dur * (W - PAD.r - PAD.l);
887+ ctx.strokeStyle = strikeColor(which); ctx.lineWidth = 2; ctx.beginPath(); ctx.moveTo(x, PAD.t); ctx.lineTo(x, H - PAD.b); ctx.stroke();
888+ }
889+ }
890+
891+ // ---------------- wiring ----------------
892+ $("playA").addEventListener("click", () => play(0, null));
893+ $("playB").addEventListener("click", () => play(1, null));
894+ $("playAB").addEventListener("click", () => play(0, 1));
895+ $("stop").addEventListener("click", stopAudio);
896+ $("slow").addEventListener("click", () => (smOn ? stopSlow() : startSlow()));
897+ for (const [id, other] of [["lisOne", "lisSum"], ["lisSum", "lisOne"]]) $(id).addEventListener("click", () => {
898+ if ($(id).getAttribute("aria-pressed") === "true") return;
899+ $(id).setAttribute("aria-pressed", "true"); $(other).setAttribute("aria-pressed", "false"); readParams(); run();
900+ });
901+ let timer = 0;
902+ for (const id of ["lx", "ly", "E0", "Kb", "Hd", "N", "sx", "sy", "rad", "dl", "tlo", "thi", "f1", "dur"]) $(id).addEventListener("input", () => { readParams(); clearTimeout(timer); timer = setTimeout(run, 400); });
903+ document.querySelectorAll("[data-preset]").forEach((b) => b.addEventListener("click", () => { setControls(PRESETS[b.dataset.preset]); run(); }));
904+ let rt = 0;
905+ window.addEventListener("resize", () => { clearTimeout(rt); rt = setTimeout(drawAll, 120); });
906+ try { window.matchMedia("(prefers-color-scheme: dark)").addEventListener("change", drawAll); new MutationObserver(drawAll).observe(document.documentElement, { attributes: true, attributeFilter: ["data-theme"] }); } catch (e) {}
907+
908+ $("lx").value = 0.62; $("ly").value = 0.7;
909+ setControls(PRESETS.gong);
910+ run();
911+})();
912+</script>