fc06d37Bead-and-spring percussion simulator: numbl-scripted setup, Verlet/Yoshida solver with exact damping, linear modes, browser appJeremy Magland 1// Web Worker: runs the numbl script, computes linear modes, and simulates.
3import { relDiff } from "./sim/analysis";
4import { computeModes, isPlanar } from "./sim/modes";
5import { numblEig, runScript } from "./sim/runScript";
6import { simulate } from "./sim/solver";
7import type { BeadSystem, SimResult } from "./sim/types";
9const helperSources = import.meta.glob("../matlab/helpers/*.m", { query: "?raw", import: "default", eager: true }) as Record<string, string>;
10const helpers = Object.entries(helperSources).map(([path, source]) => ({ name: path.split("/").pop()!, source }));
12const MAX_MODE_BEADS = 1500;
14export interface RunPayload {
15 sys: {
16 n: number;
17 pos0: Float64Array;
18 pinned: Uint8Array;
19 si: Int32Array;
20 sj: Int32Array;
21 nDashpots: number;
22 listenBeads: Int32Array;
23 hitBeads: Int32Array;
24 method: string;
25 duration: number;
26 };
27 result: SimResult;
28 modeFreqs: Float64Array | null;
29 notes: string[];
30 check: { dt: number; windows: { t: number; rel: number }[] } | null;
31 // the same run with every hit scaled by 1e-3, its output scaled back up
32 linear: { audio: Float32Array; windows: { t: number; rel: number }[] } | null;
33}
35export interface RunRequest {
36 type: "run";
37 script: string;
38 dtCheck: boolean;
39 linear: boolean;
40}
42export type FromWorker =
43 | { type: "log"; text: string }
44 | { type: "stage"; text: string }
45 | { type: "progress"; t: number; duration: number; label: string }
46 | { type: "result"; payload: RunPayload }
47 | { type: "error"; message: string };
49const post = (m: FromWorker, transfer: Transferable[] = []) => (self as unknown as Worker).postMessage(m, transfer);
51const windowsOf = (duration: number) =>
52 [0.01, 0.03, 0.1, 0.3, 1, duration].filter((t, i, a) => t <= duration && a.indexOf(t) === i);
54self.onmessage = (ev: MessageEvent<RunRequest>) => {
55 const { script, dtCheck, linear } = ev.data;
56 try {
57 post({ type: "stage", text: "Running script" });
58 const { sys, warnings } = runScript(script, helpers, text => post({ type: "log", text }));
59 const notes = [...warnings];
61 let modes = null;
62 const nFree = sys.pinned.reduce((a, p) => a + (p ? 0 : 1), 0);
63 if (sys.modes) {
64 if (!isPlanar(sys)) notes.push("Linear modes skipped: the beads do not rest in a plane z = const.");
65 else if (nFree > MAX_MODE_BEADS) notes.push(`Linear modes skipped: ${nFree} free beads (limit ${MAX_MODE_BEADS}).`);
66 else {
67 post({ type: "stage", text: `Computing linear modes (${nFree} free beads)` });
68 modes = computeModes(sys, numblEig);
69 }
70 }
72 post({ type: "stage", text: "Simulating" });
73 const result = simulate(sys, {
74 modes,
75 onProgress: t => post({ type: "progress", t, duration: sys.duration, label: "Simulating" }),
76 });
78 let checkOut: RunPayload["check"] = null;
79 if (dtCheck && !result.aborted) {
80 post({ type: "stage", text: "Simulating again with dt/2" });
81 const half: BeadSystem = { ...sys, dt: result.dt / 2 };
82 const r2 = simulate(half, {
83 maxAnimBytes: 0,
84 onProgress: t => post({ type: "progress", t, duration: sys.duration, label: "Simulating with dt/2" }),
85 });
86 const fs = result.sampleRate;
87 const windows = windowsOf(sys.duration).map(t => ({ t, rel: relDiff(r2.audio, result.audio, Math.round(t * fs)) }));
88 checkOut = { dt: r2.dt, windows };
89 }
91 let linearOut: RunPayload["linear"] = null;
92 if (linear && !result.aborted) {
93 post({ type: "stage", text: "Simulating the small-strike reference" });
94 const s = 1e-3;
95 const small: BeadSystem = {
96 ...sys,
97 dt: result.dt,
98 hits: sys.hits.map(h => ({
99 ...h,
100 velocity: [h.velocity[0] * s, h.velocity[1] * s, h.velocity[2] * s],
101 force: [h.force[0] * s, h.force[1] * s, h.force[2] * s],
102 })),
103 };
104 const r3 = simulate(small, {
105 maxAnimBytes: 0,
106 onProgress: t => post({ type: "progress", t, duration: sys.duration, label: "Small-strike reference" }),
107 });
108 const audio = r3.audio.map(v => v / s);
109 const fs = result.sampleRate;
110 const windows = windowsOf(sys.duration).map(t => ({ t, rel: relDiff(audio, result.audio, Math.round(t * fs)) }));
111 linearOut = { audio, windows };
112 }
114 const hitBeads = Int32Array.from(new Set(sys.hits.flatMap(h => Array.from(h.beads))));
115 const payload: RunPayload = {
116 sys: {
117 n: sys.n,
118 pos0: sys.pos0,
119 pinned: sys.pinned,
120 si: sys.si,
121 sj: sys.sj,
122 nDashpots: sys.di.length,
123 listenBeads: sys.listenBeads,
124 hitBeads,
125 method: sys.method,
126 duration: sys.duration,
127 },
128 result,
129 modeFreqs: modes ? modes.freqs : null,
130 notes,
131 check: checkOut,
132 linear: linearOut,
133 };
134 const transfer: Transferable[] = [result.audio.buffer, ...result.animFrames.map(f => f.buffer)];
135 if (result.modeEnergy) transfer.push(result.modeEnergy.buffer);
136 if (linearOut) transfer.push(linearOut.audio.buffer);
137 post({ type: "result", payload }, transfer);
138 } catch (e) {
139 const err = e as { message?: string; line?: number; file?: string };
140 const where = err.line ? ` (line ${err.line})` : "";
141 post({ type: "error", message: `${err.message ?? String(e)}${where}` });
142 }
143};