ae4a407cymplay: design and play a cymsim cymbal in the browserJeremy Magland 1// Per-mode listener gains: the Rayleigh front sum and the diffracted back sum over the elements,
2// with each element's travel delay folded into a phase at the mode frequency. The pressure of
3// mode m is then Re(G_m * complex acceleration), and the unbaffled estimate is front minus back.
4//
5// This is the only part of the precompute that depends on the listener, and it is cheap enough
6// (one pass over elements times modes) to redo whenever the listener or the frequencies change,
7// which is what makes the listener and the material live controls in the design panel. The rim
8// points of the diffraction path are the mesh's own boundary nodes rather than a separate polar
9// sampling, so that the same function serves the first computation and every later one.
10import type { Listener } from "../config/schema.ts";
11import type { ModesData } from "../modes.ts";
13export interface Gains {
14 front: Float32Array; // (Nm, 2) re, im
15 back: Float32Array;
16}
18interface Geometry {
19 area: Float64Array; // (E,)
20 rFront: Float64Array; // (E,)
21 rBack: Float64Array; // (E,)
22}
24/** Element areas and the two path lengths to the listener; depends on the listener position only
25 * through the paths, so it is recomputed whenever the listener moves. */
26export function radiationGeometry(data: ModesData, listener: Listener): Geometry {
27 const { xy, tri, z0, nElements, nNodes } = data;
28 const [Lx, Ly, Lz] = listener.position;
29 const area = new Float64Array(nElements);
30 const cx = new Float64Array(nElements), cy = new Float64Array(nElements), cz = new Float64Array(nElements);
31 for (let e = 0; e < nElements; e++) {
32 const a = tri[3 * e], b = tri[3 * e + 1], c = tri[3 * e + 2];
33 const x0 = xy[2 * a], y0 = xy[2 * a + 1], x1 = xy[2 * b], y1 = xy[2 * b + 1], x2 = xy[2 * c], y2 = xy[2 * c + 1];
34 area[e] = Math.abs(0.5 * ((x1 - x0) * (y2 - y0) - (x2 - x0) * (y1 - y0)));
35 cx[e] = (x0 + x1 + x2) / 3;
36 cy[e] = (y0 + y1 + y2) / 3;
37 cz[e] = (z0[a] + z0[b] + z0[c]) / 3;
38 }
39 // rim points: the nodes of the outer ring, with their own heights
40 let R = 0;
41 for (let n = 0; n < nNodes; n++) R = Math.max(R, Math.hypot(xy[2 * n], xy[2 * n + 1]));
42 const rimX: number[] = [], rimY: number[] = [], rimZ: number[] = [], rimD: number[] = [];
43 for (let n = 0; n < nNodes; n++) {
44 const x = xy[2 * n], y = xy[2 * n + 1];
45 if (Math.hypot(x, y) < R * (1 - 1e-6)) continue;
46 rimX.push(x);
47 rimY.push(y);
48 rimZ.push(z0[n]);
49 rimD.push(Math.hypot(x - Lx, y - Ly, z0[n] - Lz));
50 }
51 const nRim = rimX.length;
52 const rFront = new Float64Array(nElements), rBack = new Float64Array(nElements);
53 for (let e = 0; e < nElements; e++) {
54 rFront[e] = Math.hypot(Lx - cx[e], Ly - cy[e], Lz - cz[e]);
55 let best = Infinity;
56 for (let j = 0; j < nRim; j++) {
57 const d = Math.hypot(cx[e] - rimX[j], cy[e] - rimY[j], cz[e] - rimZ[j]) + rimD[j];
58 if (d < best) best = d;
59 }
60 rBack[e] = best;
61 }
62 return { area, rFront, rBack };
63}
65/** Gains for the given frequencies. ``phiFactor`` rescales the stored mode shapes, which are
66 * mass-normalised at the reference density; a density change scales them by 1/sqrt(rho/rho0). */
67export function computeGains(data: ModesData, listener: Listener, freq: Float64Array, phiFactor: number,
68 geom?: Geometry, onProgress?: (frac: number) => void): Gains {
69 const g = geom ?? radiationGeometry(data, listener);
70 const { tri, phiQ, phiScale, nElements, nModes } = data;
71 const front = new Float32Array(2 * nModes), back = new Float32Array(2 * nModes);
72 const fr = new Float64Array(2 * nModes), bk = new Float64Array(2 * nModes);
73 const omega = new Float64Array(nModes);
74 for (let m = 0; m < nModes; m++) omega[m] = 2 * Math.PI * freq[m];
75 const c = listener.c_air, rhoAir = listener.rho_air, gBack = listener.back_diffraction_gain;
76 const phiE = new Float64Array(nModes);
77 for (let e = 0; e < nElements; e++) {
78 const i0 = tri[3 * e] * nModes, i1 = tri[3 * e + 1] * nModes, i2 = tri[3 * e + 2] * nModes;
79 const s = phiFactor / 3;
80 for (let m = 0; m < nModes; m++) phiE[m] = s * phiScale[m] * (phiQ[i0 + m] + phiQ[i1 + m] + phiQ[i2 + m]);
81 const cf = (rhoAir * g.area[e]) / (2 * Math.PI * g.rFront[e]);
82 const cb = (gBack * rhoAir * g.area[e]) / (2 * Math.PI * g.rBack[e]);
83 const tf = g.rFront[e] / c, tb = g.rBack[e] / c;
84 for (let m = 0; m < nModes; m++) {
85 const w = omega[m], p = phiE[m];
86 const af = w * tf, ab = w * tb;
87 const pf = cf * p, pb = cb * p;
88 fr[2 * m] += pf * Math.cos(af);
89 fr[2 * m + 1] -= pf * Math.sin(af);
90 bk[2 * m] += pb * Math.cos(ab);
91 bk[2 * m + 1] -= pb * Math.sin(ab);
92 }
93 if (onProgress && (e & 255) === 0) onProgress(e / nElements);
94 }
95 front.set(fr);
96 back.set(bk);
97 return { front, back };
98}