2 * The instrument author: writes complete instrument files for the three v1
3 * instruments (drum, dulcimer, chime) from one parameterized builder.
4 *
5 * The physics and its numbers are soundbox's, taken from the same equations
6 * the numbl pipeline compiles there (the test suite holds the two projects
7 * together; see scripts/test-node.ts). What this builder adds is the form:
8 * every grid size and parameter reaches the WGSL as a pipeline override
9 * constant, so one file covers every setting of every parameter, and the
10 * whole precompute-vs-tuning split is carried by the parameter classes in
11 * the file rather than by host code that knows what a drum is.
12 *
13 * The precompute is the reciprocal room recording: run the air alone,
14 * kicked at the microphone, and after every step project the field onto
15 * each mode's coupling pattern (a product of interpolated sines under the
16 * object's radiation profile; for a string, one more pattern for its
17 * bridge patch). Discrete reciprocity makes each projection row that
18 * pattern's impulse response at the microphone. The render solves every
19 * object mode in closed form (the leapfrog restricted to a sine eigenmode
20 * is a two-term scalar recursion) and filters each kernel with a complex
21 * one-pole per root, one GPU thread per (root, kernel) pair, tiled over
22 * time with a reduction per tile.
23 *
24 * The recorded mode set is geometric: sine modes up to `mode_margin` times
25 * what the air grid can spatially resolve across the object
26 * (m <= mode_margin * L / (2h)). Nothing about the strike or the tuning
27 * enters, so the bake depends on geometry and accuracy parameters alone.
28 */
29import type { InstrumentFile, MarkerSpec, ParamSpec, PassSpec, VisualSpec } from '../src/format.ts';
31const WG = 64;
33/** Thread index linearized over a possibly two-dimensional workgroup grid
34 * (the x dimension caps at 65535 workgroups). */
35const IDX = 'let i = gid.x + gid.y * 4194240u;';
37const MAIN = `@compute @workgroup_size(${WG})
38fn main(@builtin(global_invocation_id) gid: vec3<u32>) {
39 ${IDX}`;
41const ov = (names: [string, 'u32' | 'f32'][]): string =>
42 names.map(([n, t]) => `override ${n}: ${t};`).join('\n');
44const buf = (bindings: [string, 'read' | 'read_write'][]): string =>
45 bindings
46 .map(
47 ([n, m], i) =>
48 `@group(0) @binding(${i}) var<storage, ${m}> ${n === '__tick' ? 'tick' : n === '__events' ? 'ev' : n}: array<${n === '__tick' ? 'u32' : 'f32'}>;`,
49 )
50 .join('\n');
52export type ObjectKind = 'membrane' | 'string' | 'bar';
53export type BodyKind = 'shell' | 'box' | 'none';
55export interface BuildOptions {
56 name: string;
57 description: string;
58 object: ObjectKind;
59 body: BodyKind;
60 /** Geometry parameter defaults overriding the built-in ones. */
61 defaults?: Record<string, number>;
62}
64/* ------------------------------------------------------------- params --- */
66const P = (
67 key: string,
68 label: string,
69 unit: string,
70 value: number,
71 min: number,
72 max: number,
73 cls: 'geometry' | 'accuracy' | 'tuning',
74 hint?: string,
75): ParamSpec => ({ key, label, unit, value, min, max, class: cls, hint });
77const bodyParams: Record<BodyKind, ParamSpec[]> = {
78 none: [],
79 shell: [
80 P('shellr', 'shell radius', 'm', 0.17, 0.08, 0.3, 'geometry', 'A hair larger than the head over it.'),
81 P('shelld', 'shell depth', 'm', 0.3, 0.05, 0.4, 'geometry', 'The cavity under the head; deeper is a lower air resonance.'),
82 P('thick', 'wall thickness', 'm', 0.016, 0.008, 0.05, 'geometry'),
83 P('closed', 'closed bottom', '1', 1, 0, 1, 'geometry', '1 puts a plate across the bottom; 0 leaves it open.'),
84 P('absorb', 'surface absorption', '1/s', 800, 0, 4000, 'geometry', 'A lossy skin over the shell surfaces; 0 lets the cavity ring hard.'),
85 ],
86 box: [
87 P('boxl', 'box length', 'm', 0.66, 0.3, 0.9, 'geometry'),
88 P('boxw', 'box width', 'm', 0.16, 0.08, 0.4, 'geometry'),
89 P('boxd', 'box depth', 'm', 0.08, 0.03, 0.25, 'geometry', 'How far the box hangs below its top plate at z = 0.'),
90 P('thick', 'wall thickness', 'm', 0.016, 0.008, 0.05, 'geometry'),
91 P('holer', 'sound hole radius', 'm', 0.035, 0, 0.12, 'geometry', '0 seals the box; the hole is what lets the cavity speak.'),
92 P('holex', 'sound hole x', 'm', -0.05, -0.3, 0.3, 'geometry'),
93 P('absorb', 'surface absorption', '1/s', 800, 0, 4000, 'geometry'),
94 ],
95};
97const objectParams: Record<ObjectKind, { geometry: ParamSpec[]; tuning: ParamSpec[] }> = {
98 membrane: {
99 geometry: [
100 P('L1', 'head width', 'm', 0.32, 0.15, 0.5, 'geometry', 'Extent along x.'),
101 P('L2', 'head depth', 'm', 0.32, 0.15, 0.5, 'geometry', 'Extent along y.'),
102 P('oz', 'head height', 'm', 0.02, 0, 0.15, 'geometry', 'Height of the head plane above z = 0.'),
103 ],
104 tuning: [
105 P('f0', 'pitch', 'Hz', 110, 60, 300, 'tuning', 'The (1,1) fundamental; the tension follows from it and the size.'),
106 P('t60', 'decay t60', 's', 1.5, 0.2, 6, 'tuning', 'Seconds to decay 60 dB.'),
107 P('sig1', 'brightness decay', 'm^2/s', 0.0005, 0, 0.01, 'tuning', 'Frequency-dependent damping; the slap of the attack fades first.'),
108 P('strikex', 'strike x', '1', 0.35, 0.1, 0.9, 'tuning', 'Where the stick lands along x, as a fraction.'),
109 P('strikey', 'strike y', '1', 0.5, 0.1, 0.9, 'tuning'),
110 P('striker', 'stick width', 'm', 0.04, 0.01, 0.1, 'tuning', 'Radius of the struck patch.'),
111 P('vel', 'strike speed', 'm/s', 1.5, 0.1, 5, 'tuning', 'Linear equations: sets scale only.'),
112 P('gair', 'radiation gain', '1', 1, 0, 2, 'tuning'),
113 ],
114 },
115 string: {
116 geometry: [
117 P('L1', 'speaking length', 'm', 0.6, 0.3, 0.9, 'geometry'),
118 P('oy', 'across', 'm', 0, -0.15, 0.15, 'geometry'),
119 P('oz', 'height', 'm', 0.02, 0.005, 0.15, 'geometry', 'Height of the string above z = 0.'),
120 ],
121 tuning: [
122 P('f0', 'pitch', 'Hz', 294, 80, 880, 'tuning', 'The fundamental; the wave speed follows, cs = 2 L f0.'),
123 P('B', 'inharmonicity', '1', 0.0001, 0, 0.002, 'tuning', 'Bending stiffness as the inharmonicity coefficient.'),
124 P('t60', 'decay t60', 's', 4, 0.2, 8, 'tuning'),
125 P('sig1', 'brightness decay', 'm^2/s', 0.005, 0, 0.05, 'tuning'),
126 P('pluckpos', 'pluck position', '1', 0.22, 0.08, 0.92, 'tuning'),
127 P('amp', 'pluck height', 'm', 0.002, 0.0002, 0.005, 'tuning'),
128 P('gair', 'line radiation gain', '1', 0.3, 0, 2, 'tuning'),
129 P('gbridge', 'bridge drive gain', '1', 1, 0, 2, 'tuning', 'The pull at the far end driving a patch above the plate.'),
130 ],
131 },
132 bar: {
133 geometry: [
134 P('L1', 'length', 'm', 0.3, 0.15, 0.5, 'geometry'),
135 P('oy', 'across', 'm', 0, -0.15, 0.15, 'geometry'),
136 P('oz', 'height', 'm', 0.05, 0.005, 0.15, 'geometry'),
137 ],
138 tuning: [
139 P('f1', 'pitch', 'Hz', 780, 200, 2000, 'tuning', 'The first partial; higher ones land at n^2 f1.'),
140 P('t60', 'decay t60', 's', 6, 0.5, 12, 'tuning'),
141 P('sig1', 'brightness decay', 'm^2/s', 0.001, 0, 0.02, 'tuning'),
142 P('strikepos', 'strike position', '1', 0.3, 0.1, 0.9, 'tuning'),
143 P('striker', 'mallet width', 'm', 0.025, 0.01, 0.08, 'tuning'),
144 P('vel', 'strike speed', 'm/s', 1, 0.1, 5, 'tuning'),
145 P('gair', 'radiation gain', '1', 1, 0, 2, 'tuning'),
146 ],
147 },
148};
150/* ------------------------------------------------------------ derived --- */
152function derivedOf(object: ObjectKind): Record<string, string> {
153 const air: Record<string, string> = {
154 c0: '343',
155 cfl: '0.5',
156 ny: 'nx / 2',
157 nz: 'nx / 2',
158 npts: 'nx * ny * nz',
159 h: 'Lx / nx',
160 Ly: 'ny * h',
161 Lz: 'nz * h',
162 dt: 'cfl * h / (sqrt(3) * c0)',
163 ir_steps: 'max(1, round(ir_seconds / dt))',
164 mic_index:
165 'clamp(round((mic_x + Lx/2) / h - 0.5), 0, nx - 1)' +
166 ' + nx * (clamp(round((mic_y + Ly/2) / h - 0.5), 0, ny - 1)' +
167 ' + ny * clamp(round((mic_z + Lz/2) / h - 0.5), 0, nz - 1))',
168 };
169 if (object === 'membrane') {
170 return {
171 ...air,
172 // The grid, sized for the worst case any tuning slider can reach, so
173 // no tuning move ever invalidates it (soundbox's makeRectGrid).
174 c_worst: '2 * f0__max / sqrt(1/(L1*L1) + 1/(L2*L2))',
175 h_min: 'sqrt(2 * (c_worst*c_worst*dt*dt + 4*sig1__max*dt)) / 0.95',
176 n1: 'max(5, min(256, floor(L1 / h_min) + 1))',
177 n2: 'max(5, min(256, floor(L2 / h_min) + 1))',
178 h1: 'L1 / (n1 - 1)',
179 h2: 'L2 / (n2 - 1)',
180 // The geometric mode set: what the air grid resolves, times a margin.
181 mx: 'min(n1 - 2, ceil(mode_margin * L1 / (2 * h)))',
182 my: 'min(n2 - 2, ceil(mode_margin * L2 / (2 * h)))',
183 n_rows: 'mx * my',
184 };
185 }
186 const stiff: Record<string, string> =
187 object === 'string'
188 ? {
189 cs_worst: '2 * L1 * f0__max',
190 kap2_worst: 'B__max * (cs_worst * cs_worst * L1 * L1) / (pi * pi)',
191 }
192 : {
193 kap_worst: '2 * f1__max * (L1 * L1) / pi',
194 kap2_worst: 'kap_worst * kap_worst',
195 cs_worst: '0',
196 };
197 return {
198 ...air,
199 ...stiff,
200 a_worst: 'cs_worst*cs_worst*dt*dt + 4*sig1__max*dt',
201 h_min: 'sqrt((a_worst + sqrt(a_worst*a_worst + 16*kap2_worst*dt*dt)) / 2) / 0.95',
202 n1: 'max(5, min(1024, floor(L1 / h_min) + 1))',
203 h1: 'L1 / (n1 - 1)',
204 mx: 'n1 - 2', // line objects have few modes; record all of them
205 n_rows: object === 'string' ? 'mx + 1' : 'mx', // + the bridge patch row
206 };
207}
209/* ---------------------------------------------------- medium (bodies) --- */
211const bodyWgsl: Record<BodyKind, { code: string; consts: [string, 'f32'][] }> = {
212 none: { code: ' let shell = 0.0;', consts: [] },
213 shell: {
214 code: ` // The shell: a cylindrical ring from z = -shelld to 0, optional bottom.
215 let r = sqrt(x*x + y*y);
216 let cyl_o = 0.5 * (1.0 - tanh((r - shellr) / sw));
217 let cyl_i = 0.5 * (1.0 - tanh((r - (shellr - thick)) / sw));
218 let cyl_z = 0.5 * (1.0 - tanh((abs(z + shelld/2.0) - shelld/2.0) / sw));
219 let bot = 0.5 * (1.0 - tanh((abs(z + shelld - thick/2.0) - thick/2.0) / sw));
220 let shell = min(1.0, max(0.0, cyl_o - cyl_i) * cyl_z + closed * (cyl_i * bot));`,
221 consts: [['shellr', 'f32'], ['shelld', 'f32'], ['thick', 'f32'], ['closed', 'f32']],
222 },
223 box: {
224 code: ` // The box: the outer solid minus its cavity, top plate at z = 0, with a
225 // round sound hole cut out of the top plate only.
226 let box_ox = 0.5 * (1.0 - tanh((abs(x) - boxl/2.0) / sw));
227 let box_oy = 0.5 * (1.0 - tanh((abs(y) - boxw/2.0) / sw));
228 let box_oz = 0.5 * (1.0 - tanh((abs(z + boxd/2.0) - boxd/2.0) / sw));
229 let box_ix = 0.5 * (1.0 - tanh((abs(x) - (boxl/2.0 - thick)) / sw));
230 let box_iy = 0.5 * (1.0 - tanh((abs(y) - (boxw/2.0 - thick)) / sw));
231 let box_iz = 0.5 * (1.0 - tanh((abs(z + boxd/2.0) - (boxd/2.0 - thick)) / sw));
232 var shell = max(0.0, box_ox * box_oy * box_oz - box_ix * box_iy * box_iz);
233 let hole_r = sqrt((x - holex)*(x - holex) + y*y);
234 let hole = 0.5 * (1.0 - tanh((hole_r - holer) / sw));
235 let top = 0.5 * (1.0 - tanh((abs(z + thick/2.0) - thick/2.0) / sw));
236 shell = max(0.0, shell - hole * top);`,
237 consts: [
238 ['boxl', 'f32'], ['boxw', 'f32'], ['boxd', 'f32'], ['thick', 'f32'],
239 ['holer', 'f32'], ['holex', 'f32'],
240 ],
241 },
242};
244function mediumShader(object: ObjectKind, body: BodyKind): { code: string; constants: string[] } {
245 const b = bodyWgsl[body];
246 const absorb = body === 'none' ? '' : ' + absorb * (4.0 * shell * (1.0 - shell))';
247 const profConsts: [string, 'f32'][] =
248 object === 'membrane'
249 ? [['oz', 'f32'], ['L1', 'f32'], ['L2', 'f32']]
250 : [['oy', 'f32'], ['oz', 'f32'], ['L1', 'f32']];
251 const prof =
252 object === 'membrane'
253 ? ` // The head radiates from its plane: a thin sheet at z = oz over its
254 // extent, cross-section normalized to one.
255 let aw = 2.0 * h;
256 let sk = exp(-((z - oz)*(z - oz)) / (aw*aw)) / (sqrt(PI) * aw);
257 let exx = 0.5 * (1.0 + tanh((L1/2.0 - abs(x)) / (1.5*h)));
258 let exy = 0.5 * (1.0 + tanh((L2/2.0 - abs(y)) / (1.5*h)));
259 prof[i] = sk * exx * exy * w;`
260 : ` // The object radiates from its line: a tube two cells wide around
261 // (y = oy, z = oz), fading at its ends, cross-section normalized to one.
262 let aw = 2.0 * h;
263 let tube = exp(-((y - oy)*(y - oy) + (z - oz)*(z - oz)) / (aw*aw)) / (PI * aw * aw);
264 let span = 0.5 * (1.0 + tanh((L1/2.0 - abs(x)) / (1.5*h)));
265 prof[i] = tube * span * w;`;
266 const bprof =
267 object === 'string'
268 ? ` // ...and its bridge force drives a patch of air just above the plate,
269 // around the x = +L/2 end (a stand-in for a moving top plate).
270 let pw = 1.5 * h;
271 let pr = exp(-((x - L1/2.0)*(x - L1/2.0) + (y - oy)*(y - oy)) / (0.04*0.04)) / (PI * 0.04 * 0.04);
272 let skin = exp(-(z - pw)*(z - pw) / (pw*pw)) / (sqrt(PI) * pw);
273 bprof[i] = pr * skin * w;`
274 : '';
275 const buffers: [string, 'read_write'][] = [
276 ['wall', 'read_write'], ['sig', 'read_write'], ['prof', 'read_write'],
277 ];
278 if (object === 'string') buffers.push(['bprof', 'read_write']);
279 const consts: [string, 'u32' | 'f32'][] = [
280 ['npts', 'u32'], ['nx', 'u32'], ['ny', 'u32'],
281 ['h', 'f32'], ['Lx', 'f32'], ['Ly', 'f32'], ['Lz', 'f32'],
282 ...b.consts,
283 ...(body === 'none' ? [] : ([['absorb', 'f32']] as [string, 'f32'][])),
284 ...profConsts,
285 ];
286 const code = `// The medium: the body's rigid mask, the absorption, and where the object
287// couples into the air. One thread per air cell.
288const PI = 3.14159265358979;
289${ov(consts)}
290${buf(buffers)}
291${MAIN}
292 if (i >= npts) { return; }
293 let ix = i % nx;
294 let iy = (i / nx) % ny;
295 let iz = i / (nx * ny);
296 let x = -Lx/2.0 + (f32(ix) + 0.5) * h;
297 let y = -Ly/2.0 + (f32(iy) + 0.5) * h;
298 let z = -Lz/2.0 + (f32(iz) + 0.5) * h;
299 let sw = 0.45 * h;
300${b.code}
301 let w = 1.0 - shell;
302 wall[i] = w;
303 // Open boundary at the domain edge (a quadratic absorbing ramp over the
304 // outer 0.15 Ly of each face), plus the body's lossy surface skin.
305 let spw = 0.15 * Ly;
306 let dx = max(0.0, spw - (Lx/2.0 - abs(x))) / spw;
307 let dy = max(0.0, spw - (Ly/2.0 - abs(y))) / spw;
308 let dz = max(0.0, spw - (Lz/2.0 - abs(z))) / spw;
309 let d = max(max(dx, dy), dz);
310 sig[i] = 9000.0 * d * d${absorb};
311${prof}
312${bprof}
313}
314`;
315 return { code, constants: consts.map(([n]) => n) };
316}
318/* ------------------------------------------- the rest of the precompute --- */
320/** Interpolated sine tables: the pattern the rad coupling op applies is the
321 * node sine linearly interpolated to each air column, so the projection
322 * must use exactly that, not the continuous sine. */
323function sinesShader(axis: 'x' | 'y'): { code: string; constants: string[] } {
324 const [n, m, L, hobj, nair, Lair] =
325 axis === 'x' ? ['n1', 'mx', 'L1', 'h1', 'nx', 'Lx'] : ['n2', 'my', 'L2', 'h2', 'ny', 'Ly'];
326 const consts: [string, 'u32' | 'f32'][] = [
327 [n, 'u32'], [m, 'u32'], [nair, 'u32'],
328 [L, 'f32'], [hobj, 'f32'], [Lair, 'f32'], ['h', 'f32'],
329 ];
330 const code = `// Sine mode ${axis === 'x' ? 'm' : 'n'} of the object grid, interpolated to air
331// ${axis} centers: s${axis}[(mode-1) * ${nair} + i${axis}]. One thread per table entry.
332const PI = 3.14159265358979;
333${ov(consts)}
334${buf([['table', 'read_write']])}
335${MAIN}
336 if (i >= ${m} * ${nair}) { return; }
337 let mi = i % ${m};
338 let ia = i / ${m};
339 let xw = -${Lair}/2.0 + (f32(ia) + 0.5) * h;
340 let u = (xw + ${L}/2.0) / ${hobj};
341 let j = i32(floor(u));
342 var v = 0.0;
343 if (j >= 0 && j < i32(${n}) - 1) {
344 let f = u - f32(j);
345 let mm = f32(mi + 1u);
346 let s0 = sin(mm * PI * f32(j) / f32(${n} - 1u));
347 let s1 = sin(mm * PI * f32(j + 1) / f32(${n} - 1u));
348 v = (1.0 - f) * s0 + f * s1;
349 }
350 table[mi * ${nair} + ia] = v;
351}
352`;
353 return { code, constants: consts.map(([q]) => q) };
354}
356const impulseShader = {
357 code: `// The reciprocal kick: one unit impulse at the microphone cell, entering
358// the way a source enters the first step. One thread.
359${ov([['mic_index', 'u32'], ['dt', 'f32']])}
360${buf([['p', 'read_write'], ['wall', 'read'], ['sig', 'read']])}
361@compute @workgroup_size(1)
362fn main() {
363 p[mic_index] = dt * dt * wall[mic_index] / (1.0 + sig[mic_index] * dt);
364}
365`,
366 constants: ['mic_index', 'dt'],
367};
369const stepShader = {
370 code: `// One leapfrog step of the air: the wall-masked Laplacian (each stencil
371// face scaled by the mask at the neighbour it reads, the rigid condition in
372// divergence form) fused with the update. One thread per air cell.
373${ov([['npts', 'u32'], ['nx', 'u32'], ['ny', 'u32'], ['nz', 'u32'], ['h', 'f32'], ['dt', 'f32'], ['c0', 'f32']])}
374${buf([['p', 'read'], ['pm', 'read'], ['pn', 'read_write'], ['wall', 'read'], ['sig', 'read']])}
375fn face(ix: i32, iy: i32, iz: i32, pc: f32) -> f32 {
376 if (ix < 0 || ix >= i32(nx) || iy < 0 || iy >= i32(ny) || iz < 0 || iz >= i32(nz)) {
377 return -pc; // outside the domain: p = 0, open
378 }
379 let k = u32(ix) + nx * (u32(iy) + ny * u32(iz));
380 return wall[k] * (p[k] - pc);
381}
382${MAIN}
383 if (i >= npts) { return; }
384 let ix = i32(i % nx);
385 let iy = i32((i / nx) % ny);
386 let iz = i32(i / (nx * ny));
387 let pc = p[i];
388 let s = face(ix - 1, iy, iz, pc) + face(ix + 1, iy, iz, pc)
389 + face(ix, iy - 1, iz, pc) + face(ix, iy + 1, iz, pc)
390 + face(ix, iy, iz - 1, pc) + face(ix, iy, iz + 1, pc);
391 let lp = wall[i] * s / (h * h);
392 let sd = sig[i] * dt;
393 pn[i] = wall[i] * ((2.0 * pc - (1.0 - sd) * pm[i] + (c0 * dt) * (c0 * dt) * lp) / (1.0 + sd));
394}
395`,
396 constants: ['npts', 'nx', 'ny', 'nz', 'h', 'dt', 'c0'],
397};
399function planeShader(profName: string, outName: string): { code: string; constants: string[] } {
400 return {
401 code: `// Contract the field against the ${profName === 'bprof' ? 'bridge patch' : 'radiation'} profile over z:
402// one thread per air column.
403${ov([['nx', 'u32'], ['ny', 'u32'], ['nz', 'u32']])}
404${buf([[outName, 'read_write'], ['p', 'read'], [profName, 'read']])}
405${MAIN}
406 if (i >= nx * ny) { return; }
407 var s = 0.0;
408 for (var iz = 0u; iz < nz; iz++) {
409 let k = i + nx * ny * iz;
410 s += ${profName}[k] * p[k];
411 }
412 ${outName}[i] = s;
413}
414`,
415 constants: ['nx', 'ny', 'nz'],
416 };
417}
419const xcShader = {
420 code: `// Contract the plane against every x sine row: tmp[m + mx * iy].
421${ov([['nx', 'u32'], ['ny', 'u32'], ['mx', 'u32']])}
422${buf([['tmp', 'read_write'], ['plane', 'read'], ['sx', 'read']])}
423${MAIN}
424 if (i >= mx * ny) { return; }
425 let m = i % mx;
426 let iy = i / mx;
427 var s = 0.0;
428 for (var ix = 0u; ix < nx; ix++) {
429 s += sx[m * nx + ix] * plane[ix + nx * iy];
430 }
431 tmp[i] = s;
432}
433`,
434 constants: ['nx', 'ny', 'mx'],
435};
437/** The last contraction writes this step's sample of every kernel row. */
438function projShader(object: ObjectKind): { code: string; constants: string[] } {
439 if (object === 'membrane') {
440 return {
441 code: `// Finish each pattern's projection against its y sine row and write this
442// step's sample: kernels[q * ir_steps + t], q = (m-1) + mx*(n-1).
443${ov([['ny', 'u32'], ['mx', 'u32'], ['my', 'u32'], ['ir_steps', 'u32']])}
444${buf([['kernels', 'read_write'], ['tmp', 'read'], ['sy', 'read'], ['__tick', 'read']])}
445${MAIN}
446 if (i >= mx * my) { return; }
447 let m = i % mx;
448 let n = i / mx;
449 var s = 0.0;
450 for (var iy = 0u; iy < ny; iy++) {
451 s += sy[n * ny + iy] * tmp[m + mx * iy];
452 }
453 let t = tick[0];
454 if (t < ir_steps) {
455 kernels[i * ir_steps + t] = s;
456 }
457}
458`,
459 constants: ['ny', 'mx', 'my', 'ir_steps'],
460 };
461 }
462 return {
463 code: `// A line object's y factor is one everywhere: finish by summing tmp over
464// iy and write this step's sample of every kernel row.
465${ov([['ny', 'u32'], ['mx', 'u32'], ['ir_steps', 'u32']])}
466${buf([['kernels', 'read_write'], ['tmp', 'read'], ['__tick', 'read']])}
467${MAIN}
468 if (i >= mx) { return; }
469 var s = 0.0;
470 for (var iy = 0u; iy < ny; iy++) {
471 s += tmp[i + mx * iy];
472 }
473 let t = tick[0];
474 if (t < ir_steps) {
475 kernels[i * ir_steps + t] = s;
476 }
477}
478`,
479 constants: ['ny', 'mx', 'ir_steps'],
480 };
481}
483const bridgeProjShader = {
484 code: `// The bridge patch is a single fixed pattern: sum its plane and write the
485// sample into the last kernel row (row mx). One thread.
486${ov([['nx', 'u32'], ['ny', 'u32'], ['mx', 'u32'], ['ir_steps', 'u32']])}
487${buf([['kernels', 'read_write'], ['bplane', 'read'], ['__tick', 'read']])}
488@compute @workgroup_size(1)
489fn main() {
490 var s = 0.0;
491 for (var k = 0u; k < nx * ny; k++) {
492 s += bplane[k];
493 }
494 let t = tick[0];
495 if (t < ir_steps) {
496 kernels[mx * ir_steps + t] = s;
497 }
498}
499`,
500 constants: ['nx', 'ny', 'mx', 'ir_steps'],
501};
503/* ------------------------------------------------------------- render --- */
505/** Shared complex and slot helpers for the slots shaders. A slot is one
506 * characteristic root ready to filter one kernel row:
507 * [cRe, cIm, zRe, zIm, row, start], six f32 each. */
508const SLOT_HELPERS = `
509fn cmul(ar: f32, ai: f32, br: f32, bi: f32) -> vec2f {
510 return vec2f(ar * br - ai * bi, ar * bi + ai * br);
511}
512// c = a * gain * (z - 2 + 1/z) = a * gain * (z-1)^2 / z, the acceleration
513// series, with z-1 supplied stably (the naive z - 1 cancels).
514fn acc_scale(p: vec4f, zm1Re: f32, gain: f32) -> vec4f {
515 let m2 = p.z * p.z + p.w * p.w;
516 if (m2 == 0.0) { return vec4f(0.0); }
517 let sq = cmul(zm1Re, p.w, zm1Re, p.w);
518 let num = cmul(sq.x, sq.y, p.z, -p.w);
519 let c = cmul(p.x, p.y, gain * num.x / m2, gain * num.y / m2);
520 return vec4f(c.x, c.y, p.z, p.w);
521}
522// c = a * gain * z, the next-displacement series (what a string's bridge
523// slope is read from).
524fn next_scale(p: vec4f, gain: f32) -> vec4f {
525 let c = cmul(p.x, p.y, gain * p.z, gain * p.w);
526 return vec4f(c.x, c.y, p.z, p.w);
527}
528`;
530const WRITE_SLOT = `
531fn write_slot(idx: u32, p: vec4f, row: f32, start: f32) {
532 slots[idx * 6u] = p.x;
533 slots[idx * 6u + 1u] = p.y;
534 slots[idx * 6u + 2u] = p.z;
535 slots[idx * 6u + 3u] = p.w;
536 slots[idx * 6u + 4u] = row;
537 slots[idx * 6u + 5u] = start;
538}
539`;
541/**
542 * The mode's recursion solved in closed form, written to avoid the f32
543 * cancellations the naive discriminant carries. Restricted to eigenvalue
544 * lam, the leapfrog is a[k+1] = cA a[k] + cB a[k-1] with
545 * cA = (2 - x) / (1 + d), x = dt^2 w2 + 2 sig1 dt lam,
546 * cB = (-(1 - d) + 2 sig1 dt lam) / (1 + d), d = (6.91 / t60) dt,
547 * and the quantities the roots need reduce to
548 * zRe = (1 - x/2)/(1+d), |z|^2 = (1 - g)/(1+d), g = d + 2 sig1 dt lam,
549 * zIm^2 (1+d)^2 = N = dt^2 w2 - g d - x^2/4,
550 * every term small and positive-definite where it matters, so the mode's
551 * angle keeps full relative f32 precision (the naive cA^2 + 4 cB loses it).
552 * N > 0 is the underdamped (audible) case; N <= 0 gives two real roots.
553 * Excitation: a strike is a velocity (a[0] = 0, a[-1] = b), a pluck a held
554 * displacement (a[0] = a[-1] = b).
555 *
556 * Leaves in scope: p0, p1 (displacement-series roots as [aRe, aIm, zRe,
557 * zIm], zeroed when absent) and zm10, zm11 (each root's z - 1 real part,
558 * stably).
559 */
560const rootsBlock = (excite: 'strike' | 'pluck'): string => `
561 let d = (6.91 / t60) * dt;
562 let opd = 1.0 + d;
563 let g = d + 2.0 * sig1 * dt * lam;
564 let x = dt * dt * w2 + 2.0 * sig1 * dt * lam;
565 let zRe = (1.0 - 0.5 * x) / opd;
566 let zm1Re = -(0.5 * x + d) / opd;
567 let nn = dt * dt * w2 - g * d - 0.25 * x * x;
568 var p0 = vec4f(0.0);
569 var p1 = vec4f(0.0);
570 var zm10 = 0.0;
571 var zm11 = 0.0;
572 if (nn > 0.0) {
573 let zIm = sqrt(nn) / opd;
574 let mag2 = (1.0 - g) / opd;
575${excite === 'strike'
576 ? ` let aRe = 0.0;
577 let aIm = b * mag2 / zIm;`
578 : ` let aRe = b;
579 let aIm = b * ((0.5 * x - g) / opd) / zIm;`}
580 p0 = vec4f(aRe, aIm, zRe, zIm);
581 zm10 = zm1Re;
582 } else {
583 let sq = sqrt(-nn) / opd;
584 let r1 = zRe + sq;
585 let r2 = zRe - sq;
586 if (r1 != r2 && r1 != 0.0 && r2 != 0.0) {
587 let det = 1.0 / r1 - 1.0 / r2;
588${excite === 'strike'
589 ? ` let a1 = b / det;
590 let a2 = -a1;`
591 : ` let a1 = (b - b / r2) / det;
592 let a2 = b - a1;`}
593 p0 = vec4f(a1, 0.0, r1, 0.0);
594 p1 = vec4f(a2, 0.0, r2, 0.0);
595 zm10 = zm1Re + sq;
596 zm11 = zm1Re - sq;
597 }
598 }`;
600/** Event-row accessors for the tuning parameters a shader needs. */
601const evAccess = (tuning: ParamSpec[], needed: string[]): string =>
602 needed
603 .map((k) => {
604 const j = tuning.findIndex((p) => p.key === k);
605 if (j < 0) throw new Error(`event param '${k}' not in tuning list`);
606 return ` let ${k} = ev[base + ${1 + j}u];`;
607 })
608 .join('\n');
610/** The per-axis excitation sums: for a strike, sum_i sin(m pi i/(n-1)) *
611 * exp(-((i h - pos)/r)^2); a pluck replaces the Gaussian by the triangle.
612 * One thread per (event, mode index). */
613function excitationShader(
614 object: ObjectKind,
615 tuning: ParamSpec[],
616 axis: 'x' | 'y',
617): { code: string; constants: string[] } {
618 const [n, m, L, hobj] = axis === 'x' ? ['n1', 'mx', 'L1', 'h1'] : ['n2', 'my', 'L2', 'h2'];
619 let access: string;
620 let term: string;
621 if (object === 'string') {
622 access = evAccess(tuning, ['pluckpos', 'amp']);
623 term = ` let xi = f32(k) * ${hobj};
624 let ta = xi / (pluckpos * ${L});
625 let tb = (${L} - xi) / (${L} - pluckpos * ${L});
626 let u0 = amp * 0.5 * (ta + tb - abs(ta - tb));`;
627 } else {
628 const pos = object === 'bar' ? 'strikepos' : axis === 'x' ? 'strikex' : 'strikey';
629 access = evAccess(tuning, [pos, 'striker']);
630 term = ` let xi = f32(k) * ${hobj};
631 let dx = (xi - ${pos} * ${L}) / striker;
632 let u0 = exp(-dx * dx);`;
633 }
634 const consts: [string, 'u32' | 'f32'][] = [
635 [n, 'u32'], [m, 'u32'], ['ev_stride', 'u32'], [hobj, 'f32'], [L, 'f32'],
636 ];
637 const code = `// Project this axis's excitation profile onto each recorded sine mode:
638// b${axis}[e * ${m} + (mode-1)]. One thread per (event, mode).
639const PI = 3.14159265358979;
640${ov(consts)}
641${buf([['bout', 'read_write'], ['__events', 'read']])}
642${MAIN}
643 if (i >= ${m} * arrayLength(&ev) / ev_stride) { return; }
644 let mi = i % ${m};
645 let e = i / ${m};
646 let base = e * ev_stride;
647${access}
648 let mm = f32(mi + 1u);
649 var s = 0.0;
650 for (var k = 1u; k < ${n} - 1u; k++) {
651${term}
652 s += sin(mm * PI * f32(k) / f32(${n} - 1u)) * u0;
653 }
654 bout[i] = s;
655}
656`;
657 return { code, constants: consts.map(([q]) => q) };
658}
660/** The slots shader: every recorded mode of every event solved in closed
661 * form and scaled into filter slots. */
662function slotsShader(object: ObjectKind, tuning: ParamSpec[]): { code: string; constants: string[] } {
663 if (object === 'membrane') {
664 const consts: [string, 'u32' | 'f32'][] = [
665 ['mx', 'u32'], ['my', 'u32'], ['n_rows', 'u32'], ['n1', 'u32'], ['n2', 'u32'],
666 ['ev_stride', 'u32'],
667 ['h1', 'f32'], ['h2', 'f32'], ['L1', 'f32'], ['L2', 'f32'], ['dt', 'f32'],
668 ];
669 const code = `// One thread per (event, recorded mode): the mode's excitation amplitude,
670// its recursion roots, and the acceleration-series slots that drive its
671// kernel row.
672const PI = 3.14159265358979;
673${ov(consts)}
674${buf([['slots', 'read_write'], ['__events', 'read'], ['bx', 'read'], ['by', 'read']])}
675${SLOT_HELPERS}
676${WRITE_SLOT}
677${MAIN}
678 let n_events = arrayLength(&ev) / ev_stride;
679 if (i >= n_events * n_rows) { return; }
680 let e = i / n_rows;
681 let q = i % n_rows;
682 let base = e * ev_stride;
683 let start = ev[base];
684${evAccess(tuning, ['f0', 't60', 'sig1', 'vel', 'gair'])}
685 let mi = q % mx;
686 let ni = q / mx;
687 let m = f32(mi + 1u);
688 let n = f32(ni + 1u);
689 let b = bx[e * mx + mi] * by[e * my + ni] * (-(dt * vel)) * 4.0 / (f32(n1 - 1u) * f32(n2 - 1u));
690 let s1 = sin(m * PI / (2.0 * f32(n1 - 1u)));
691 let s2 = sin(n * PI / (2.0 * f32(n2 - 1u)));
692 let lam = (4.0 / (h1 * h1)) * s1 * s1 + (4.0 / (h2 * h2)) * s2 * s2;
693 let cm = 2.0 * f0 / sqrt(1.0 / (L1 * L1) + 1.0 / (L2 * L2));
694 let w2 = cm * cm * lam;
695${rootsBlock('strike')}
696 let ga = gair / (dt * dt);
697 write_slot(i * 2u, acc_scale(p0, zm10, ga), f32(q), start);
698 write_slot(i * 2u + 1u, acc_scale(p1, zm11, ga), f32(q), start);
699}
700`;
701 return { code, constants: consts.map(([q]) => q) };
702 }
703 if (object === 'bar') {
704 const consts: [string, 'u32' | 'f32'][] = [
705 ['mx', 'u32'], ['n1', 'u32'], ['ev_stride', 'u32'],
706 ['h1', 'f32'], ['L1', 'f32'], ['dt', 'f32'],
707 ];
708 const code = `// One thread per (event, mode) of the bar: the tension-free stiff scheme,
709// partials at n^2 f1 up to the grid's own dispersion.
710const PI = 3.14159265358979;
711${ov(consts)}
712${buf([['slots', 'read_write'], ['__events', 'read'], ['bx', 'read']])}
713${SLOT_HELPERS}
714${WRITE_SLOT}
715${MAIN}
716 let n_events = arrayLength(&ev) / ev_stride;
717 if (i >= n_events * mx) { return; }
718 let e = i / mx;
719 let mi = i % mx;
720 let base = e * ev_stride;
721 let start = ev[base];
722${evAccess(tuning, ['f1', 't60', 'sig1', 'vel', 'gair'])}
723 let m = f32(mi + 1u);
724 let b = bx[i] * (-(dt * vel)) * 2.0 / f32(n1 - 1u);
725 let s1 = sin(m * PI / (2.0 * f32(n1 - 1u)));
726 let lam = (4.0 / (h1 * h1)) * s1 * s1;
727 let kap = 2.0 * f1 * (L1 * L1) / PI;
728 let w2 = kap * kap * lam * lam;
729${rootsBlock('strike')}
730 let ga = gair / (dt * dt);
731 write_slot(i * 2u, acc_scale(p0, zm10, ga), f32(mi), start);
732 write_slot(i * 2u + 1u, acc_scale(p1, zm11, ga), f32(mi), start);
733}
734`;
735 return { code, constants: consts.map(([q]) => q) };
736 }
737 const consts: [string, 'u32' | 'f32'][] = [
738 ['mx', 'u32'], ['n1', 'u32'], ['ev_stride', 'u32'],
739 ['h1', 'f32'], ['L1', 'f32'], ['dt', 'f32'],
740 ];
741 const code = `// One thread per (event, mode) of the string: the stiff damped string's
742// roots, its radiation slots, and the bridge slots its arriving slope
743// drives (the last kernel row, row mx).
744const PI = 3.14159265358979;
745${ov(consts)}
746${buf([['slots', 'read_write'], ['__events', 'read'], ['bx', 'read']])}
747${SLOT_HELPERS}
748${WRITE_SLOT}
749${MAIN}
750 let n_events = arrayLength(&ev) / ev_stride;
751 if (i >= n_events * mx) { return; }
752 let e = i / mx;
753 let mi = i % mx;
754 let base = e * ev_stride;
755 let start = ev[base];
756${evAccess(tuning, ['f0', 'B', 't60', 'sig1', 'gair', 'gbridge'])}
757 let m = f32(mi + 1u);
758 let b = bx[i] * 2.0 / f32(n1 - 1u);
759 let s1 = sin(m * PI / (2.0 * f32(n1 - 1u)));
760 let lam = (4.0 / (h1 * h1)) * s1 * s1;
761 let cs = 2.0 * L1 * f0;
762 let kap2 = B * (cs * L1 / PI) * (cs * L1 / PI);
763 let w2 = cs * cs * lam + kap2 * lam * lam;
764${rootsBlock('pluck')}
765 let ga = gair / (dt * dt);
766 write_slot(i * 4u, acc_scale(p0, zm10, ga), f32(mi), start);
767 write_slot(i * 4u + 1u, acc_scale(p1, zm11, ga), f32(mi), start);
768 // The bridge: the far end is pinned, so the arriving slope is
769 // -u[n-2]/h1, and the generated step's term is
770 // -(gbridge cs^2 / L) bprof fb; together, +gbridge cs^2/(L h1) times
771 // sin(m pi (n1-2)/(n1-1)) times the next-displacement series.
772 let gb = (gbridge * cs * cs / (L1 * h1)) * sin(m * PI * f32(n1 - 2u) / f32(n1 - 1u));
773 write_slot(i * 4u + 2u, next_scale(p0, gb), f32(mx), start);
774 write_slot(i * 4u + 3u, next_scale(p1, gb), f32(mx), start);
775}
776`;
777 return { code, constants: consts.map(([q]) => q) };
778}
780const synthShader = {
781 code: `// One thread per slot: a complex one-pole over this tile of samples, the
782// kernel row as input while it lasts and silence after (the closed-form
783// tail). State persists across tiles.
784${ov([['tile', 'u32'], ['total_steps', 'u32'], ['ir_steps', 'u32'], ['n_slots', 'u32']])}
785${buf([['slots', 'read'], ['state', 'read_write'], ['kernels', 'read'], ['partial', 'read_write'], ['__tick', 'read']])}
786${MAIN}
787 if (i >= n_slots) { return; }
788 let cRe = slots[i * 6u];
789 let cIm = slots[i * 6u + 1u];
790 let zRe = slots[i * 6u + 2u];
791 let zIm = slots[i * 6u + 3u];
792 let row = u32(slots[i * 6u + 4u]);
793 let start = u32(slots[i * 6u + 5u]);
794 var wRe = state[i * 2u];
795 var wIm = state[i * 2u + 1u];
796 let t0 = tick[0] * tile;
797 for (var s = 0u; s < tile; s++) {
798 let t = t0 + s;
799 if (t >= total_steps) { break; }
800 var hv = 0.0;
801 if (t >= start && t - start < ir_steps) {
802 hv = kernels[row * ir_steps + (t - start)];
803 }
804 let nRe = zRe * wRe - zIm * wIm + hv;
805 let nIm = zRe * wIm + zIm * wRe;
806 wRe = nRe;
807 wIm = nIm;
808 partial[i * tile + s] = cRe * wRe - cIm * wIm;
809 }
810 state[i * 2u] = wRe;
811 state[i * 2u + 1u] = wIm;
812}
813`,
814 constants: ['tile', 'total_steps', 'ir_steps', 'n_slots'],
815};
817const reduceShader = {
818 code: `// One thread per sample of the tile: sum every slot's contribution.
819${ov([['tile', 'u32'], ['total_steps', 'u32'], ['n_slots', 'u32']])}
820${buf([['partial', 'read'], ['trace', 'read_write'], ['__tick', 'read']])}
821${MAIN}
822 if (i >= tile) { return; }
823 let t = tick[0] * tile + i;
824 if (t >= total_steps) { return; }
825 var s = 0.0;
826 for (var q = 0u; q < n_slots; q++) {
827 s += partial[q * tile + i];
828 }
829 trace[t] = s;
830}
831`,
832 constants: ['tile', 'total_steps', 'n_slots'],
833};
835/* ----------------------------------------------------------- assembly --- */
837export function buildInstrument(opts: BuildOptions): InstrumentFile {
838 const { object, body } = opts;
839 const op = objectParams[object];
840 const params: ParamSpec[] = [
841 { ...P('nx', 'air cells along x', '1', 96, 32, 256, 'accuracy',
842 'The air grid; ny = nz = nx/2. The simulation is honest to about c/(8h). Must be even.'), step: 2 },
843 P('Lx', 'domain length', 'm', 1.0, 0.6, 2.0, 'geometry',
844 'The air box around the instrument, twice as long in x as in y and z.'),
845 P('ir_seconds', 'room response length', 's', 0.5, 0.1, 2.0, 'accuracy',
846 'How long the recorded room responses run; a room still ringing past this is truncated.'),
847 ...(object === 'membrane'
848 ? [P('mode_margin', 'mode set margin', '1', 2.0, 1.0, 4.0, 'accuracy',
849 'Head modes are recorded up to this multiple of what the air grid spatially resolves.')]
850 : []),
851 ...op.geometry,
852 ...bodyParams[body],
853 P('mic_x', 'microphone x', 'm', 0.15, -0.45, 0.45, 'geometry'),
854 P('mic_y', 'microphone y', 'm', 0.1, -0.2, 0.2, 'geometry'),
855 P('mic_z', 'microphone z', 'm', 0.12, -0.2, 0.2, 'geometry'),
856 ...op.tuning,
857 ];
858 for (const [k, v] of Object.entries(opts.defaults ?? {})) {
859 const p = params.find((q) => q.key === k);
860 if (!p) throw new Error(`default for unknown parameter '${k}'`);
861 p.value = v;
862 }
863 const tuning = op.tuning;
864 const evStride = 1 + tuning.length;
866 const shaders: Record<string, { code: string; constants: string[] }> = {
867 medium: mediumShader(object, body),
868 sines_x: sinesShader('x'),
869 impulse: impulseShader,
870 step: stepShader,
871 plane: planeShader('prof', 'plane'),
872 xc: xcShader,
873 proj: projShader(object),
874 excite_x: excitationShader(object, tuning, 'x'),
875 slots: slotsShader(object, tuning),
876 synth: synthShader,
877 reduce: reduceShader,
878 };
879 if (object === 'membrane') {
880 shaders.sines_y = sinesShader('y');
881 shaders.excite_y = excitationShader(object, tuning, 'y');
882 }
883 if (object === 'string') {
884 shaders.bplane = planeShader('bprof', 'bplane');
885 shaders.bridge_proj = bridgeProjShader;
886 }
888 const preBuffers = [
889 { name: 'wall', count: 'npts' },
890 { name: 'sig', count: 'npts' },
891 { name: 'prof', count: 'npts' },
892 ...(object === 'string' ? [{ name: 'bprof', count: 'npts' }] : []),
893 { name: 'sx', count: 'mx * nx' },
894 ...(object === 'membrane' ? [{ name: 'sy', count: 'my * ny' }] : []),
895 { name: 'p', count: 'npts' },
896 { name: 'pm', count: 'npts' },
897 { name: 'pn', count: 'npts' },
898 { name: 'plane', count: 'nx * ny' },
899 ...(object === 'string' ? [{ name: 'bplane', count: 'nx * ny' }] : []),
900 { name: 'tmp', count: 'mx * ny' },
901 { name: 'kernels', count: 'n_rows * ir_steps' },
902 ];
903 const mediumOut = ['wall', 'sig', 'prof', ...(object === 'string' ? ['bprof'] : [])];
904 const preSetup: PassSpec[] = [
905 { shader: 'medium', buffers: mediumOut, dispatch: 'npts' },
906 { shader: 'sines_x', buffers: ['sx'], dispatch: 'mx * nx' },
907 ...(object === 'membrane'
908 ? [{ shader: 'sines_y', buffers: ['sy'], dispatch: 'my * ny' }]
909 : []),
910 { shader: 'impulse', buffers: ['p', 'wall', 'sig'], dispatch: '1' },
911 ];
912 const preLoop: PassSpec[] = [
913 { shader: 'plane', buffers: ['plane', 'p', 'prof'], dispatch: 'nx * ny' },
914 { shader: 'xc', buffers: ['tmp', 'plane', 'sx'], dispatch: 'mx * ny' },
915 object === 'membrane'
916 ? { shader: 'proj', buffers: ['kernels', 'tmp', 'sy', '__tick'], dispatch: 'n_rows' }
917 : { shader: 'proj', buffers: ['kernels', 'tmp', '__tick'], dispatch: 'mx' },
918 ...(object === 'string'
919 ? [
920 { shader: 'bplane', buffers: ['bplane', 'p', 'bprof'], dispatch: 'nx * ny' },
921 { shader: 'bridge_proj', buffers: ['kernels', 'bplane', '__tick'], dispatch: '1' },
922 ]
923 : []),
924 { shader: 'step', buffers: ['p', 'pm', 'pn', 'wall', 'sig'], dispatch: 'npts' },
925 ];
927 const slotsPer = object === 'string' ? 4 : 2;
928 const slotsCount = object === 'membrane' ? 'n_events * n_rows * 2' : `n_events * mx * ${slotsPer}`;
929 const renBuffers = [
930 { name: 'kernels', count: 'n_rows * ir_steps' },
931 { name: 'bx', count: 'n_events * mx' },
932 ...(object === 'membrane' ? [{ name: 'by', count: 'n_events * my' }] : []),
933 { name: 'slots', count: 'n_slots * 6' },
934 { name: 'state', count: 'n_slots * 2' },
935 { name: 'partial', count: 'n_slots * tile' },
936 { name: 'trace', count: 'total_steps' },
937 ];
938 const MIC: [number, number, number] = [0.3, 0.85, 1.0];
939 const OUTLINE: [number, number, number] = [1.0, 0.85, 0.3];
940 const markers: MarkerSpec[] = [
941 object === 'membrane'
942 ? {
943 kind: 'rect',
944 corner: ['-L1/2', '-L2/2', 'oz'],
945 du: ['L1', '0', '0'],
946 dv: ['0', 'L2', '0'],
947 color: OUTLINE,
948 label: 'head',
949 }
950 : {
951 kind: 'line',
952 from: ['-L1/2', 'oy', 'oz'],
953 to: ['L1/2', 'oy', 'oz'],
954 color: OUTLINE,
955 label: object,
956 },
957 { kind: 'point', at: ['mic_x', 'mic_y', 'mic_z'], color: MIC, label: 'mic' },
958 ];
959 const visual: VisualSpec = {
960 dims: ['nx', 'ny', 'nz'],
961 extents: ['Lx', 'Ly', 'Lz'],
962 buffers: [
963 { name: 'wall', count: 'npts' },
964 { name: 'sig', count: 'npts' },
965 { name: 'prof', count: 'npts' },
966 ...(object === 'string' ? [{ name: 'bprof', count: 'npts' }] : []),
967 ],
968 setup: [{ shader: 'medium', buffers: mediumOut, dispatch: 'npts' }],
969 volumes: [
970 // The body is the medium shader's own mask, so what the viewer shows
971 // is what the solver sees, grid effects included.
972 { buffer: 'wall', style: 'solid', invert: true, color: [0.62, 0.48, 0.34] },
973 { buffer: 'prof', style: 'cloud', color: [1.0, 0.45, 0.2] },
974 ...(object === 'string'
975 ? [{ buffer: 'bprof', style: 'cloud' as const, color: [0.35, 0.65, 1.0] as [number, number, number] }]
976 : []),
977 ],
978 markers,
979 partials: { buffer: 'slots', stride: 6, zRe: 2, zIm: 3 },
980 };
982 const renSetup: PassSpec[] = [
983 { shader: 'excite_x', buffers: ['bx', '__events'], dispatch: 'n_events * mx' },
984 ...(object === 'membrane'
985 ? [{ shader: 'excite_y', buffers: ['by', '__events'], dispatch: 'n_events * my' }]
986 : []),
987 object === 'membrane'
988 ? { shader: 'slots', buffers: ['slots', '__events', 'bx', 'by'], dispatch: 'n_events * n_rows' }
989 : { shader: 'slots', buffers: ['slots', '__events', 'bx'], dispatch: 'n_events * mx' },
990 ];
992 return {
993 instrawave: 1,
994 name: opts.name,
995 description: opts.description,
996 params,
997 derived: derivedOf(object),
998 shaders,
999 precompute: {
1000 buffers: preBuffers,
1001 setup: preSetup,
1002 loop: { count: 'ir_steps', batch: 512, rotate: [['p', 'pm', 'pn']], body: preLoop },
1003 outputs: ['kernels'],
1004 },
1005 render: {
1006 buffers: renBuffers,
1007 inputs: ['kernels'],
1008 eventParams: tuning.map((p) => p.key),
1009 derived: { ev_stride: `${evStride}`, n_slots: slotsCount },
1010 setup: renSetup,
1011 loop: {
1012 count: 'n_tiles',
1013 batch: 8,
1014 body: [
1015 { shader: 'synth', buffers: ['slots', 'state', 'kernels', 'partial', '__tick'], dispatch: 'n_slots' },
1016 { shader: 'reduce', buffers: ['partial', 'trace', '__tick'], dispatch: 'tile' },
1017 ],
1018 },
1019 trace: 'trace',
1020 },
1021 visual,
1022 };
1023}