/ proof-of-concept / steelsim41
Sign in
proof-of-concept / steelsim41
steelsim41 / THEORY.md
107 lines · 27.3 KBPreviewCodeBlameHistoryRaw
1# Model and numerical method
3This note describes the model behind the simulator: a thin metal plate of arbitrary outline and spatially varying thickness, with a stress-free curved shape (the bow and the bell of a cymbal), either hanging freely or sitting on a stand, struck once with a drumstick. The plate rings like a cymbal or a gong. Its vibration amplitude, several thicknesses, is beyond the range where a linear plate model is accurate, so we solve the geometrically nonlinear plate equations of Foppl and von Karman with a thickness field, keep the in-plane inertia, discretize on a regular grid, and step explicitly on the GPU. Section 1 gives the continuous model, Section 2 the support and the rigid-body motion, Section 3 the strike, Sections 4 and 5 the discretization in space and time, Section 6 the radiation model, and Section 7 the limitations. The method is that of steelsim31, the hanging steel sheet, extended to a thickness field, a general outline, and a mount; the parts that carry over are repeated here so that this note stands alone.
5## 1. Plate model
7The plate occupies a region $\Omega$ of the plane, given by the design script as a bounding box $[-L_x/2, L_x/2] \times [-L_y/2, L_y/2]$ and an outline function (a rectangle, a circle or an ellipse with or without a center hole, or any other test). It has a thickness $h(x,y)$, Young's modulus $E$, Poisson ratio $\nu$, and density $\rho$; the material constants are uniform, the thickness is not. Its motion is described by the in-plane displacements $u(x,y,t)$ and $v(x,y,t)$ and the transverse displacement $w(x,y,t)$, all measured from a stress-free initial shape $w_0(x,y)$, the profile of the design (below). The membrane strains in the Marguerre form are
9$$
10\varepsilon_{xx} = u_x + \tfrac12\left[(w_0 + w)_x^2 - w_{0,x}^2\right], \qquad
11\varepsilon_{yy} = v_y + \tfrac12\left[(w_0 + w)_y^2 - w_{0,y}^2\right],
12$$
14$$
15\gamma_{xy} = u_y + v_x + (w_0 + w)_x (w_0 + w)_y - w_{0,x} w_{0,y},
16$$
18so that $u = v = w = 0$ carries no stress. The membrane forces per unit length are $N_{xx} = C(\varepsilon_{xx} + \nu\varepsilon_{yy})$, $N_{yy} = C(\varepsilon_{yy} + \nu\varepsilon_{xx})$, and $N_{xy} = G\gamma_{xy}$, with $C = Eh/(1-\nu^2)$ and $G = Eh/(2(1+\nu))$, both now functions of position through $h$. Writing $W = w_0 + w$ for the total shape and $D = Eh^3/(12(1-\nu^2))$ for the local bending stiffness, the equations of motion are
20$$
21\rho h\, u_{tt} = \partial_x N_{xx} + \partial_y N_{xy}, \qquad
22\rho h\, v_{tt} = \partial_x N_{xy} + \partial_y N_{yy},
23$$
25$$
26\rho h\, w_{tt} = -\mathcal{B}[w] + \partial_x(N_{xx} W_x + N_{xy} W_y) + \partial_y(N_{xy} W_x + N_{yy} W_y) - 2\sigma_0\rho h\, w_t + 2\sigma_1 \rho h\, \nabla^2 w_t ,
27$$
29where $\mathcal{B}[w]$ is the variable-stiffness bending operator, the Euler-Lagrange derivative of the bending energy
31$$
32U_b = \frac12 \int_\Omega D(x,y)\left[w_{xx}^2 + w_{yy}^2 + 2\nu\, w_{xx}w_{yy} + 2(1-\nu)\, w_{xy}^2\right] dA ,
33$$
35which reduces to $D\nabla^4 w$ where $D$ is constant. Bending acts on the deviation $w$ from the initial shape; the membrane forces act through the total slope. For $w_0 = 0$ and constant $h$ these are the von Karman equations. The quadratic terms in the strains couple bending to stretching and produce the cymbal-like behavior: the pitch glide of a strongly excited plate and the transfer of energy from the low modes toward high frequencies, which is what turns the initial tone of a struck cymbal into a wash. The in-plane inertia terms $\rho h u_{tt}$ and $\rho h v_{tt}$ are usually dropped in favor of an Airy stress function, since in-plane waves are much faster than bending waves. We keep them because they make the whole system explicit and local, which suits the GPU, at the cost of a time step set by the in-plane wave speed (Section 5).
37**Thickness field.** A real cymbal is lathed to a taper, thick at the bell and thin at the edge, and the taper is a large part of what distinguishes a crash from a ride. The thickness enters the model in three places: the mass per unit area $\rho h$, the membrane stiffnesses $C$ and $G$, which are proportional to $h$, and the bending stiffness $D$, which is proportional to $h^3$. Note that a plate of twice the thickness has twice the bending frequencies but the same in-plane frequencies; the tests check both. The thickness must be positive over the plate, and it should vary slowly on the scale of the grid.
39**Profile.** The stress-free shape $w_0$ is given by the design script as a function of position: a paraboloidal bow plus a cosine bell for the cymbal presets, a cylindrical bow or a saddle for the sheet-metal presets, or anything else. A curved plate is a shallow shell, and its behavior differs qualitatively from a flat plate: the flat plate has a symmetric, purely hardening nonlinearity, whereas a shell has an asymmetric one, may snap between two curved states, and oscillates more slowly at large amplitude as it passes through the soft flat configuration. Only the slopes of $w_0$ enter the strains, so its mean is immaterial, but its tilt is not: the von Karman strains are not invariant under rotation (Section 2), and a tilted stress-free shape would add a spurious membrane coupling linear in the tilt. The mass-weighted least-squares plane is therefore subtracted from every profile. Bending acts on the deviation $w$ only, so the profile carries no residual stress; this is the idealization of a permanent set, whereas a real hammered and lathed cymbal does carry some. The profile should be smooth on the scale of the grid, and its slopes should be small compared with one for the shallow-shell treatment to hold. This is the main point of tension between the model and a real cymbal: the bell of a crash rises 25 mm over a radius of 50 mm, a slope of 0.5 or more, and the bow itself has a slope of 0.1 to 0.2 near the edge. The presets use somewhat flatter shapes, the interface reports the largest slope and warns above 0.1, and the china preset, with its upturned edge, is the least trustworthy of them (Section 7).
41**Losses.** The two-parameter loss model is the standard one in plate sound synthesis. A bending mode of angular frequency $\omega$ and wavenumber $k$ decays at the rate $\sigma_0 + \sigma_1 k^2$, and for the flat plate $k^2 = \omega/\kappa$ with $\kappa = \sqrt{D/\rho h}$, so $\sigma_0$ sets the low-frequency decay and $\sigma_1$ makes high frequencies decay faster in proportion to frequency. The $\sigma_1$ term is implemented as a Laplacian of the velocity and does not depend on the thickness; where the plate is thicker, $\kappa$ is larger and a given frequency corresponds to a smaller $k^2$, so the frequency-dependent loss is weaker there. The in-plane motion carries an additional frequency-independent loss $\sigma_{\mathrm{in}}$. The stretching modes of the plate lie at several kilohertz, couple linearly to the transverse motion once the plate is curved, and would otherwise ring for the whole low-frequency decay time; a material loss factor of about $10^{-3}$ at 5 kHz gives a rate of this order. The three values are set by the design script (steel: $\sigma_0 = 1.15\ \mathrm{s^{-1}}$, $\sigma_1 = 5\times10^{-4}\ \mathrm{m^2/s}$, $\sigma_{\mathrm{in}} = 20\ \mathrm{s^{-1}}$ as in steelsim31; the bronze presets ring longer). None of these values are measured; they are chosen to sound plausible.
43## 2. Support and rigid-body motion
45The whole edge is free, including the edge of a center hole: the bending moment and the Kirchhoff shear vanish, and so does the in-plane traction. Rather than imposing these conditions explicitly, the discretization of Section 4 is derived from a discrete energy, and free edges are then its natural boundary conditions.
47The plate is supported in one of two ways, chosen by the design script.
49**Free plate.** A sheet hanging from a thread is treated as free, and every external load is made self-equilibrated: the reaction to the stick is distributed over the plate as the inertial load of a rigid body, which is the inertia-relief equilibrium of a free body. The plate then gains no net momentum or angular momentum, which stands in for the thread taking them up. The stick force is known at every step, so its reaction is applied directly as the rigid acceleration field that balances the force and its moments about the center (Section 3). Once the plate rings freely, no load remains that could move it as a whole, apart from rounding and the small non-invariance of the von Karman energy under finite rotation, and an optional mass-weighted projection of the six rigid-body modes (three translations, three rotations) every 32 steps removes that drift; it is off by default.
51**Mounted plate.** A cymbal sits on a stand: a rod through the center hole, with a felt washer below and usually another above, held by a wing nut. At audio frequencies the felt is a soft spring compared with the local stiffness of the plate, so it does not clamp the cymbal; it fixes the average position and tilt of the cymbal, and it lets the cymbal rock on the stand after a hit. We model the mount as a spring of total stiffness $K_m$ and a dashpot of total rate $C_m$ to ground, spread over the plate nodes within a radius $r_m$ of the mount point in proportion to their quadrature weight, and acting on all three displacement components. A uniform translation of the plate then feels exactly $K_m$, and a tilt feels a restoring moment $K_m \langle r^2 \rangle$ set by the second moment of the mount patch, so the finite radius of the felt is what holds the cymbal level; the tests check both. With the default values ($r_m = 20$ mm, $K_m = 2\times10^4$ N/m, $C_m = 5$ N s/m) a 16 inch bronze cymbal of about 1 kg bounces on the mount at roughly 20 Hz and rocks at a few hertz, and a hit near the edge starts a visible rock of several millimeters at the edge that decays over a second or so, as on a real stand. The mount takes the reaction of the stick directly, so no inertia relief is applied to a mounted plate, and the rigid-body projection is not applied either.
53A rigid rotation is not exactly a null mode of the von Karman energy: the strains contain $\tfrac12 W_x^2$, so a rotation $\theta$ about an in-plane axis produces a strain $\theta^2/2$ plus, on a curved plate, a term $\theta\, w_{0,y}$ linear in the tilt. Both are relieved by an in-plane displacement, $v = -\theta w_0$ for the linear term, which the free edges allow and which the in-plane waves establish within a fraction of a millisecond, far faster than the rock of the cymbal. The mount's in-plane spring resists this relief only weakly ($K_m$ against a membrane stiffness $Eh$ of order $10^8$ N/m), so the rocking of a mounted cymbal does not stress it spuriously. This is the reason the loads on a free plate have to be self-equilibrated, since a slow drift in rotation there would have no such relief through the mount and would build up. Gravity is dropped: the in-plane stress it induces is a strain of order $10^{-7}$.
55## 3. Excitation
57**Drumstick strike.** The stick tip is a free mass $M$ (default 30 g, the effective mass of a stick in a wrist stroke rather than its full mass) arriving with speed $v_0$ (default 3 m/s) at a point of the plate, given as fractions of the bounding box. The contact patch consists of the plate nodes within the tip radius (5 mm) plus the nearest one. Its contact coordinate is the mass-weighted mean transverse displacement $w_c$ of these nodes, the penetration is $\delta = w_c - x$ with $x$ the tip position, and the force is the Hertz law $F = K\,\delta^{\alpha}$ for $\delta > 0$ and zero otherwise, with $\alpha = 1.5$ for a rounded tip and $K = 5\times10^{8}\ \mathrm{N/m^{1.5}}$ for wood on metal. The force is spread over the patch in proportion to the node masses, which carry the local thickness. Because that distribution is exactly conjugate to the mean displacement, the contact has a potential energy and the coupled system is Hamiltonian apart from the losses. On a free plate the force is balanced over the whole plate by the inertia-relief field $a_k = F\,(1/M_s + y_k y_c/I_x + x_k x_c/I_y)$, where $M_s$ is the mass of the plate, $I_x$ and $I_y$ its moments of inertia about the in-plane axes, and $(x_c, y_c)$ the mass-weighted center of the patch, so that the total force and the two moments of the load vanish (Section 2). The field uses the force of the previous step, a lag of one microsecond whose net impulse over the contact is zero. On a mounted plate no such field is applied. The stick obeys $M\ddot x = F$ and is stepped together with the plate, so within the hit it may chatter, releasing and re-engaging as the plate runs ahead of it and comes back, which on a thin plate happens within a fraction of a millisecond. Once the contact has been released for 5 ms the stick is withdrawn and takes no further part. Without this rule a light plate on a soft mount, which the hit sets rocking as a whole, swings back after some tens of milliseconds and strikes the stick a second time; a drummer lifts the stick well before that. The recording starts at the impact, with the plate at rest in its stress-free shape. Note that a thin plate is a soft target: its driving-point impedance is $8\sqrt{D\rho h}$, about 60 N s/m for 0.8 mm steel and about the same for 1 mm bronze, so the contact lasts well under a millisecond and the peak force is a few hundred newtons. The excitation is therefore broadband, and the sound is a short bright attack followed by the ring.
59## 4. Spatial discretization
61The bounding box is covered by a regular grid of $n_x \times n_y$ nodes with spacings $h_x$ and $h_y$ (8 mm by default), and the outline selects the part of the grid that is the plate. All three displacements live at the nodes, and so does the thickness. The forces are the negative gradient of a discrete potential energy $U = U_b + U_m + U_s$ (bending, membrane, mount), so the semi-discrete system is Hamiltonian, its linear part has a symmetric stiffness matrix, and the free edge needs no special formulas.
63**Outline.** A cell belongs to the plate when its center lies inside the outline, and a node when it is a corner of such a cell; the edge is therefore a staircase of whole cells, and a hole is simply a region whose cells are absent. Every rule below refers to this membership: a second difference exists where the node and both of its neighbors belong to the plate, a cell contributes membrane energy and the mixed bending term when its four corners do, and the weight of a node is the fraction of its four cells that belong to the plate (1 inside, 1/2 on an edge, 1/4 at a corner of a full rectangle). Because the forces are the gradient of the discrete energy, the staircase edge is free, the stiffness stays symmetric, and exactly three rigid modes remain in each block, hole or no hole. Nodes outside the plate are frozen and carry no mass; a node that would be left with no cell of its own is removed. Deciding membership by cell centers rather than by nodes keeps the area unbiased; a node-based mask loses about half a cell all around the edge and raises every frequency by several percent.
65**Bending with a thickness field.** We evaluate $w_{xx}$ and $w_{yy}$ by three-point second differences at the nodes and $w_{xy}$ by the four-point difference at cell centers, and we integrate the bending energy with the node weights, the bending stiffness $D_k = E h_k^3/(12(1-\nu^2))$ at the nodes and $D_c$ at the cells, where the cell thickness is the mean of its four corners:
67$$
68U_b = \frac{h_x h_y}{2}\left[\sum_k \omega_k D_k \left(a_k^2 + b_k^2 + 2\nu a_k b_k\right) + 2(1-\nu)\sum_c D_c\, c_c^2\right],
69$$
71with $a_k$, $b_k$ the second differences at node $k$ (absent where undefined) and $c_c$ the mixed difference at cell $c$. The force at a node is the gradient of this sum, which the kernel assembles from the weighted moments $M^a_k = \omega_k D_k (a_k + \nu b_k)$, $M^b_k = \omega_k D_k (b_k + \nu a_k)$, and $M^c_c = D_c c_c$ by applying the second-difference stencils once more. Where $D$ is constant the interior stencil is the standard 13-point biharmonic; with a varying $D$ the same gather produces the discrete form of the variable-stiffness operator, and the edges are handled exactly as the energy dictates. The three rigid modes $w = 1, x, y$ lie in the null space for any thickness field.
73**Membrane.** The cells are bilinear (Q4) finite elements with 2 by 2 Gauss integration, using the strains of Section 1 with $w_0$ and $w$ interpolated bilinearly and the stiffnesses $C$ and $G$ of the cell's mean thickness. One-point integration would be cheaper but admits hourglass modes, which the nonlinear coupling would excite. On the GPU the force of each cell on its four nodes is computed once per cell within a workgroup and gathered by the nodes.
75**Mass.** The mass matrix is lumped: the mass of a node is $\rho h_k h_x h_y \omega_k$, so the total mass is the exact integral of $\rho h$ over the staircase plate.
77**Mount.** The spring of Section 2 is a diagonal term $\tfrac12 \sum_k K_k (u_k^2 + v_k^2 + w_k^2)$ with $K_k = K_m \omega_k / \sum_{m} \omega_m$ over the mount patch, and the dashpot is the same distribution of $C_m$ applied to the velocities.
79**Accuracy.** Two errors limit the accuracy of the spatial scheme. The first is the treatment of the edge. For the rectangle the scheme converges at second order: the first seven frequencies of a completely free square plate are low by about 3.5 percent at 10 cells per side, 0.8 percent at 20, and 0.3 percent at 30, compared with the reference values of Gorman. The staircase edge of a curved outline is only first-order accurate; the frequencies of a free circular plate are within about 2 percent of Leissa's values at 24 cells across the diameter and within 1 percent at 60, and the degenerate pairs, which the staircase splits, differ by 1 to 2 percent. A 16 inch cymbal at 8 mm has about 50 cells across, so its low modes are placed to within a percent or two, and the splitting of the pairs is of the same order as the splitting that the hammering of a real cymbal produces. The second error is grid dispersion, which affects every partial regardless of the outline. A plane bending wave of wavenumber $k$ along a grid axis has the discrete frequency ratio $(2 - 2\cos kh)/(kh)^2$, so a partial comes out low by 1 percent at $kh \approx 0.35$, by 3 percent at $kh \approx 0.60$, and by 10 percent at $kh \approx 1.12$. Since the bending frequency at a given $k$ is proportional to the thickness, these points move up with the thickness: for 0.8 mm steel on the 8 mm grid the 3 percent point is near 1.1 kHz and the 10 percent point near 3.6 kHz; for 1.5 mm bronze they are near 1.6 and 5.4 kHz. The interface reports the two frequencies at the mean thickness of the design. The high partials of the cymbal sound are therefore audibly flattened at 8 mm, and a finer grid is a matter of patience rather than of the method: the cost scales with the inverse cube of the spacing (Section 5).
81**Validation.** The tests check that both stiffness blocks are symmetric with exactly three rigid modes on uniform and tapered rectangles, on an ellipse, and on a disc with a hole; that the forces are the exact gradient of the discrete energy at a strongly nonlinear state, flat and bowed, uniform and tapered, free and mounted; that the loss terms and the mount's dashpot only remove energy; that doubling the thickness doubles the bending frequencies and leaves the membrane frequencies unchanged; that the mount resists a rigid translation with exactly its stiffness and a tilt with a restoring moment; and that the free-plate frequencies converge as stated. A further test runs one time step of the GPU kernel on a smooth random state of about 2.5 thicknesses amplitude on the flat and bowed rectangles, the tapered and mounted crash and china, and the twisted ellipse, and compares it with the float64 reference: the accelerations agree to about $10^{-6}$ relative, which is single-precision rounding.
83## 5. Time integration
85The system is advanced with the symplectic Euler scheme, $v^{n+1} = v^n + \Delta t\, a(x^n, v^n)$ and $x^{n+1} = x^n + \Delta t\, v^{n+1}$, in single precision on the GPU. With the velocities interpreted at half steps this is the leapfrog (Stormer-Verlet) scheme, second order and symplectic for the conservative part; the losses and the mount's dashpot are evaluated explicitly at the old velocity. The stability limit is $\Delta t < 2/\omega_{\max}$, where $\omega_{\max}$ is the largest frequency of the linearized discrete system. We estimate $\omega_{\max}$ at start-up by power iteration on a small uniform grid with the same spacing and the largest thickness of the design, which bounds the bending frequencies, and use 90 percent of the limit. The in-plane waves, with speed $c_p = \sqrt{E/\rho(1-\nu^2)}$, about 5200 m/s for steel and 3800 m/s for bronze, set this limit, and $c_p$ does not depend on the thickness: at 8 mm spacing $\Delta t \approx 1.3\ \mu$s for steel and $1.7\ \mu$s for bronze, several times smaller than bending alone would require even at the bell of a ride, and 4 s of sound take 2 to 3 million steps; at 4 mm the step halves and the node count quadruples. The step is rounded down so that an integer number of steps makes one audio sample.
87**In-plane mass scaling.** Since the in-plane waves set the step and the in-plane inertia is a small effect (the usual Airy stress function formulation drops it altogether), the interface offers a factor $s \geq 1$ by which the in-plane inertia is scaled; the default is 1. The in-plane wave speed becomes $c_p/\sqrt{s}$ and the stable step grows by up to $\sqrt{s}$, until bending sets the limit. The trade-off is that the in-plane modes of the plate, which start near 5 to 8 kHz for the presets, move down by the same factor and, once the plate is curved, are coupled to the transverse motion and to the sound. A factor of 4 halves the run time and keeps the in-plane modes above 3 kHz; it is meant for quick previews at a finer spacing rather than for the final sound.
89Single precision needs one precaution. With $w$ of a few millimeters and $\Delta t \approx 1\ \mu$s, a velocity below about $10^{-3}$ m/s produces an increment smaller than the resolution of $w$, and the rounding would act as a broadband noise source. The update of $w$ is therefore compensated: the rounding error of each addition is carried in a separate array and added back at the next step (Kahan summation). The self-test of Section 4 also checks this update.
91The state is stored twice and each step reads one copy and writes the other. A step is one compute dispatch in which each workgroup of 256 threads advances a 12 by 12 tile of nodes: the threads load the tile with a two-node halo (displacements, profile, thickness, and node weights) into shared memory, compute the membrane forces of the cells and the bending moments at the nodes once each, and then evaluate the update of the tile's nodes from shared memory; the transverse velocity, which only the loss Laplacian needs, is read from global memory to keep the workgroup within the 16 KB of shared memory that every WebGPU adapter guarantees. On grids of a few thousand nodes the cost of a step is the latency of one such workgroup, about 15 microseconds on an integrated GPU, rather than the arithmetic, so the step time barely grows until the grid exceeds several thousand nodes. A per-step uniform, addressed by a dynamic offset, carries the damping and the audio bookkeeping, so that a whole movie frame of steps is encoded in one pass, and the following frame is encoded while the GPU works on the current one. The mass-weighted sums over the plate (rigid-body components, kinetic energy, peak displacement) run in one workgroup once per frame, or every 32 steps when the optional projection of Section 2 is on.
93## 6. Radiation to the listener
95The listener model is the Rayleigh integral for a baffled plate,
97$$
98p(\mathbf r, t) = \frac{\rho_{\mathrm{air}}}{2\pi} \sum_k \frac{\Delta A_k}{r_k}\, \ddot w_k(t - r_k/c),
99$$
101where the sum runs over the nodes, $\Delta A_k$ is the node's share of the area, $r_k$ is the distance from node $k$ to the listener, and $c = 343$ m/s. The transverse velocity of every node is recorded in a ring buffer at the audio rate, and each output sample is a delayed, weighted sum over the plate, evaluated on the GPU once per movie frame for all the frame's samples, with linear interpolation between samples; the acceleration is formed as the difference between successive delayed velocities, which is the average acceleration over one sample interval. This matters in single precision: the force evaluation carries a broadband rounding error of the order of $1\ \mathrm{m/s^2}$ per node, and sampling the instantaneous acceleration would put that noise straight into the audio, whereas the averaged form suppresses it by the square root of the number of steps per sample. The common delay is discarded; only the spread of delays across the plate is kept.
103The plate is unbaffled in reality, so the sound from the two faces partly cancels at low frequencies, and the Rayleigh integral overestimates the low-frequency level; in particular the slow rock of a mounted cymbal, at a few hertz, radiates nothing audible in reality and next to nothing here either, since the acceleration involved is tiny. An accurate treatment would require solving the exterior acoustic problem, which is outside the scope of this simulator. The pressure is reported in pascals and the peak sound pressure level is displayed. Playback and the WAV file use, by default, a fixed mapping from pressure to full scale (1 Pa at 0 dBFS, shifted by a user gain in decibels), so that runs with different designs can be compared by ear; a soft clipper bounds the output and the resulting peak level is shown. Two normalized modes are also offered, to the peak or to the loudness of the sustained sound.
105## 7. Limitations
107The von Karman model assumes moderate rotations, so it is not valid for displacements that are large compared with the plate size; a strike that moves a 400 mm cymbal by a few millimeters is well within its range. The Marguerre treatment of the profile is a shallow-shell approximation, adequate for a bow of a few millimeters on a sheet and for a gentle cymbal bow, but strained by the bell of a real cymbal and by the upturned edge of a china, whose slopes are several times the range in which the approximation holds; the presets use flatter shapes than real cymbals, and the interface reports the largest slope of a design. A deep bell also behaves as a stiff cup that shifts the bell modes up and couples them to the bow in ways the shallow shell does not capture. The thickness field is treated as a slowly varying coefficient in a thin-plate model, which is appropriate for a lathed taper but not for a step or a groove on the scale of the grid. The hammering of a real cymbal, which leaves residual stresses and a slightly non-axisymmetric shape, is absent unless the script adds it (a random field through `noise2d` is one way). The loss model is empirical and its parameters are not calibrated against a real cymbal. Grid dispersion lowers the high partials by the amounts given in Section 4, and the staircase edge places the low frequencies within a few percent and splits modes that would be degenerate on a true disc. The mount is a linear spring, whereas a felt is nonlinear and the cymbal can lift off it; and the rod in the hole, which limits the in-plane motion and the rock, is only represented by the mount's in-plane spring. Air loading, which lowers the frequencies of a thin plate slightly and adds radiation damping, is not included, and the radiation model is approximate, as discussed above. Finally, single precision limits the dynamic range of the tail of the sound to roughly 100 dB below the peak.
↑↓ move↵openescclose