4ddcc1bCymbal simulator: shell-of-revolution modal model with a nonlinear energy cascade, numbl MATLAB reference, and interactive web appJeremy Magland 1// Mirror of ../matlab/cymbal_profile.m: the cymbal meridian z(r), used by the
2// 3D scene and the mode gallery. The physics reads the same geometry from the
3// MATLAB side; this is only for drawing.
4//
5// Bow (r >= bell radius): a parabola in the rim-referred coordinate, flat at
6// the rim. Bell (r < bell radius): the quartic that matches the bow's height,
7// slope and curvature at the join and comes to a horizontal apex at r = 0.
9export interface Meridian {
10 rHole: number;
11 rBell: number;
12 rRim: number;
13 z(r: number): number;
14 thickness(r: number): number;
15 apexHeight: number;
16}
18export function cymbalMeridian(
19 radius: number,
20 holeRadius: number,
21 bellRadius: number,
22 bellHeight: number,
23 bowHeight: number,
24 thicknessBell: number,
25 thicknessRim: number,
26): Meridian {
27 const a = radius;
28 const rh = Math.max(1e-4, Math.min(holeRadius, 0.5 * a));
29 const rb = Math.min(Math.max(bellRadius, 2 * rh), 0.8 * a);
30 const den = a - rb;
31 const mb = (-2 * bowHeight) / den;
32 const kb = (2 * bowHeight) / (den * den);
33 // solve the 3x3 for the quartic coefficients (B r^2 + C r^3 + E r^4)
34 const M = [
35 [rb * rb, rb ** 3, rb ** 4],
36 [2 * rb, 3 * rb * rb, 4 * rb ** 3],
37 [2, 6 * rb, 12 * rb * rb],
38 ];
39 const rhs = [-bellHeight, mb, kb];
40 const [B, C, E] = solve3(M, rhs);
41 const A0 = bowHeight + bellHeight;
42 return {
43 rHole: rh,
44 rBell: rb,
45 rRim: a,
46 apexHeight: A0,
47 z(r: number): number {
48 if (r < rb) return A0 + B * r * r + C * r ** 3 + E * r ** 4;
49 const u = (a - r) / den;
50 return bowHeight * u * u;
51 },
52 thickness(r: number): number {
53 const taper = Math.min(1, Math.max(0, (a - r) / den));
54 return thicknessRim + (thicknessBell - thicknessRim) * taper;
55 },
56 };
57}
59function solve3(M: number[][], b: number[]): [number, number, number] {
60 const A = M.map((row, i) => [...row, b[i]]);
61 for (let c = 0; c < 3; c++) {
62 let piv = c;
63 for (let r = c + 1; r < 3; r++) if (Math.abs(A[r][c]) > Math.abs(A[piv][c])) piv = r;
64 [A[c], A[piv]] = [A[piv], A[c]];
65 const d = A[c][c] || 1e-30;
66 for (let k = c; k < 4; k++) A[c][k] /= d;
67 for (let r = 0; r < 3; r++) {
68 if (r === c) continue;
69 const f = A[r][c];
70 for (let k = c; k < 4; k++) A[r][k] -= f * A[c][k];
71 }
72 }
73 return [A[0][3], A[1][3], A[2][3]];
74}