/ simulated_instruments / cymbal
Sign in
simulated_instruments / cymbal
cymbal / README.md
9.1 KBPreviewCodeBlameHistoryRaw

Cymbal#

A physically-based simulation of a struck cymbal, sampled as acoustic pressure at a point in space. The instrument is a bronze shell of revolution with a raised bell and a tapered bow; a Hertz stick contact drives its bending and stretching modes, computed by a one-dimensional finite element eigensolve per nodal-diameter order, and the geometric nonlinearity that turns a hard stroke into a crash is carried as a cascade of energy up the modal frequency ladder. The plate radiates to a listener point as a sheet of normal dipoles.

There is no microphone model: no directivity, no capsule response, no gain staging. The output is p(r, t) in pascals.

Contents#

Path What it is
matlab/ The algorithm in MATLAB syntax, run by numbl
webapp/ Interactive browser app that executes the MATLAB code client-side

The MATLAB code runs from the command line with cd matlab && numbl run cymbal_demo.m, or cymbal_validate.m to run the checks. The web app is a Vite project; see webapp/README.md.

The model#

The signal chain, with each stage in its own .m file:

meridian profile: bell, bow, thickness taper       cymbal_profile.m
  -> thin-shell elements for one Fourier order     cymbal_shell.m
      -> one eigensolve per nodal-diameter order   cymbal_basis.m
         + radiated power -> radiation damping
stroke gesture -> Hertz stick contacts             cymbal_strikes.m
  -> modal energies, and their nonlinear flow
     up the frequency ladder                       cymbal_cascade.m
      -> oscillator bank + dipole-sheet transfer   cymbal_radiate.m
          -> free-field pressure at the listener
              -> anti-aliased decimation -> PCM    cymbal_output.m

Three points are worth stating plainly, because they are where the model commits itself.

The modes come from a shell of revolution, one Fourier order at a time. A cymbal is axisymmetric, so its displacement field separates exactly in the azimuth: every mode carries n nodal diameters and varies as cos(n theta), and the two-dimensional eigenproblem breaks into one one-dimensional problem per n. Each is discretized along the meridian with two-node Love-Kirchhoff shell elements, linear in the meridional and circumferential displacements and cubic Hermite in the normal one, integrated by Gauss quadrature against the tabulated profile rather than frozen into conical frusta (cymbal_shell.m). Rotary inertia is kept, which shifts the retained frequencies by well under a percent but is what makes the mass matrix positive definite, so each order reduces through the Cholesky factor of the mass matrix and is diagonalized by svd. The curvature of the bow is what earns this formulation its keep: on a flat disc, bending and stretching are independent, and it is the meridian curvature that couples them and pulls the partials of a cymbal away from the Chladni frequencies of a plate. Shear deformation is neglected, and the mounting is taken as free-free with a point dashpot at the hole for the felt washer, which reaches the n = 0 and n = 1 modes and almost nothing else. Only the cos member of each degenerate pair is retained: that is exact for the model as posed, since every contact sits on the meridian theta = 0 where the sin member has a node. A real cymbal is not perfectly round, and the resulting slow beating of each split pair is outside this model.

The nonlinearity is carried as energy, not as a coupled integration. A cymbal struck hard is not a linear resonator. Its transverse displacement is comparable to the plate thickness, the von Karman stretching of the mid-surface couples the modes, and energy injected low by the stick migrates upward over tens of milliseconds, which is what the ear hears as a crash opening into shimmer after its attack. Integrating that directly means a dense cubic coupling tensor over hundreds of modes, so what is integrated instead is the energy of each mode (cymbal_cascade.m): each mode has an up-target a fixed frequency ratio above it and a down-target the same ratio below, and energy performs a biased random walk between them at a rate cascade_gain * om * (w/t)^2, proportional to the mode's own squared amplitude-to-thickness ratio because that is the order of the geometric nonlinearity. Steps are taken in log frequency rather than in mode index on purpose, so the swell does not depend on how densely the band happens to be resolved. Because the flux out of a mode goes as the square of its energy, a nearly empty mode passes nothing on, the front of the cascade advances rather than the whole band draining at once, and the swell stops on its own. What this buys is the instrument's central behaviour: the same cymbal pings when tapped and crashes when hit, with no velocity layers anywhere in the model. What it gives up is phase. Amplitude comes from the cascade and phase from a linear impulse response carried alongside it, so a stroke re-phases the modes it finds quiet and a coherent attack survives, but the intermodulation products of a real cascade do not. cascade_gain is the one fitted number in the model; setting it to zero leaves an exactly linear plate, which sounds like a small gong.

Radiation is a dipole sheet, computed two ways because two regimes are needed. Each surface element carries a dipole along the local normal of volume dS * ldip, where ldip is half the shorter of the flexural and acoustic wavelengths, capped at the plate radius: the distance over which the air short-circuits from one face to the other. That length is the one phenomenological quantity in the radiation model, and below the plate's coincidence frequency, about 18 kHz for 1 mm bronze, it is the flexural wavelength that sets it, which is why a cymbal is such a poor radiator in the low register and improves steadily as the cascade fills the top. For the listener, who at half a metre from a 40 cm plate is nowhere near the far field, the element dipoles are summed coherently over the azimuth with their true distances and the near-field term kept (cymbal_radiate.m). For the radiated power, which sets the radiation damping and is defined only in the far field, the azimuthal integral of the same kernel collapses onto Bessel functions and the power integral reduces to one quadrature in cos(theta) (cymbal_basis.m); the two agree to a few parts in a thousand, which the validation checks. Summing element powers incoherently, as the sibling triangle does for a rod, would not do here: for a continuous surface that estimate halves when the elements do. The multipole structure comes out of the azimuthal sum on its own. A mode with n nodal diameters is a 2^n-pole, so the high-order modes are weak until their wavelength is short, and on the symmetry axis every one of them cancels exactly.

Validation#

matlab/cymbal_validate.m (also npm run selftest in webapp/) checks the shell elements in their flat limit against the closed-form free circular plate, whose frequency equation is the vanishing of a two-by-two determinant in J_n and I_n under zero radial moment and zero Kelvin-Kirchhoff shear, together with the rigid-body count each Fourier order must have (two at n = 0, two at n = 1, none above); the mass-orthonormality and stiffness-diagonalization identities of the eigensolve; the Bessel far-field reduction against brute-force summation over a finely discretized dipole sheet; the exact cancellation of every mode with a nodal diameter on the symmetry axis; the Hertz v^(-1/5) contact-time scaling in the light-stick limit; exact energy conservation of the cascade scheme with the damping off and the top of the ladder reflecting, and its collapse onto closed-form damped exponentials when cascade_gain is zero; the far-field -6 dB per doubling of distance; exact reproducibility of the seeded stroke pattern; and pinned end-to-end diagnostics.

As in the sibling maraca, triangle, recorder and bass repositories, there is no separate Python reference implementation: the MATLAB code is the reference, and the pinned numbers were recorded from it when the model was built.

Defaults#

A 16-inch B20 bronze crash-ride: 20.3 cm radius, a 4.5 cm bell rising 1.7 cm above a bow that itself rises 2 cm from the rim, and a thickness tapering from 1.15 mm at the bell to 0.75 mm at the edge, which comes to 1.05 kg. Solving to 48 nodal diameters retains 303 modes between 49 Hz and 8.0 kHz. The default gesture is a single crash on the edge at 3 m/s with a 30 g nylon-tipped stick held 40 degrees off the surface normal: an 88 microsecond Hertz contact, a 291 N peak force, a peak plate displacement of 150 micrometres or 0.16 thicknesses, and a peak of 13.4 Pa (116 dB SPL) at the 58 cm listener before normalization. About 6% of the struck energy leaves the top of the retained band, which is the shimmer above 8 kHz the model does not carry; raising max_mode_hz and n_diameters together recovers it at the cost of the eigensolves. Setting strike_count to 8 gives a repeated figure alternating between the edge and the bow, and dropping strike_velocity to a few tenths of a metre per second turns the same instrument back into a ping, because the cascade has nothing left to work with.

↑↓ move↵openescclose