This study investigates whether residual membrane vibration measurably influences the energy transferred by a subsequent strike during djembe performance. A modal finite element model is developed in which each strike is represented by a spatially distributed contact pressure profile applied over a finite contact duration — Hertzian-inspired for the finger and slap strikes, and annular for the cupped-palm bass strike — parameterized separately for each technique. Rayleigh damping coefficients are calibrated against the exponential decay of two spectral peaks extracted from a recorded bass strike (R² ≈ 0.95).
Because the modal equations are linear while vibrational energy is quadratic in the modal coordinates, the energy of two overlapping strikes is not the sum of their isolated energies. The residual cross term quantifies their interaction. In a bass–slap–finger sequence it reaches 0.2% of peak total energy, rising to approximately 20% at strike coincidences in two 2:3 polyrhythms, in both cases exceeding the numerical superposition error by at least an order of magnitude.
Its sign is not determined by simultaneity. Two identical slap strikes applied simultaneously at mirrored positions reinforce one another, whereas a slap and a finger strike applied simultaneously at different positions partially cancel. The governing factor is whether the two strikes excite each shared mode in phase or in anti-phase, which is fixed by the strike locations and that mode’s nodal geometry. The first antisymmetric mode confirms this: in the mirrored rhythm its two strike centers straddle that mode’s nodal diameter, and the interaction it carries is suppressed by two to three orders of magnitude.
The vibration modes through which interaction occurs are likewise governed by strike geometry. The fundamental dominates when a bass strike is present, owing to its slow decay, whereas higher modes dominate when both strikes fall near the rim and couple weakly to the fundamental. Strike placement therefore determines not only the strength of interaction between successive strikes but the region of the modal spectrum in which it takes place.

Introduction

The vibration of stretched membranes constitutes one of the fundamental problems in mathematical and computational physics. The behavior of such systems is governed by the two-dimensional wave equation,​​
∂²u
∂t²
= c²∇²u,​where u(x,y,t) denotes the transverse displacement of the membrane and c =
T
ρₛ
is the wave speed determined by the membrane tension T and surface density ρₛ. Its solutions describe the propagation of transverse waves across a tensioned surface. Although analytical solutions for perfectly circular membranes have been extensively studied using Bessel functions, these solutions assume idealized conditions including uniform membrane tension, infinitesimal point excitation, negligible contact duration and the absence of damping. Real percussion instruments violate each of these assumptions, and hence require a more concrete modelling approach.​The djembe is a single-headed goblet drum renowned for a wide tonal range, capable of producing deep bass notes, resonant open tones, and bright slap sounds using only variations in hand position and striking technique. This versatility arises from the coupled interaction between the vibrating membrane, the player’s hand, and the enclosed air cavity. The behavior of a djembe membrane is further complicated by finite-area and finite-time hand contact, energy dissipation, and striking techniques that produce different spatial and temporal loading conditions.​These physical complexities make analytical solutions impractical for realistic percussion instruments. Instead, using numerical methods is an effective means of approximating the governing equations while preserving the underlying mechanics of the system. Among these, the Finite Element Method (FEM) has become one of the most widely used techniques for solving partial differential equations because it can accurately represent complex geometries, spatially varying material properties, and realistic loading conditions.​

Motivation

Personal Motivation

I have been playing hand percussion instruments, particularly the djembe and the darbuka, for approximately two years. During this time I observed that relatively small changes in strike location, hand geometry and striking technique consistently produced significantly different acoustic responses. Learning to control these subtle differences is one of the most challenging aspects of playing the instrument.
The physical origin of these differences is rather complex. This motivated an investigation into whether the distinctive sounds produced by different strikes could be explained quantitatively using computational physics rather than qualitative musical description.
Beyond understanding individual strike types, I became interested in the interaction between successive strikes during a djembe solo performance. Djembe music frequently incorporates rapid rhythmic compositions and polyrhythms in which the membrane may still be vibrating when the next strike occurs. This naturally raises the question of whether previous membrane motion influences subsequent strikes, and if so, under what temporal conditions such interactions become significant. To address this question, a computational framework is developed and applied to investigate how residual membrane vibration influences subsequent strikes during djembe performance.

Physics Motivation

From a physical perspective, membrane vibration provides an excellent example of how continuum mechanics, differential equations and numerical methods combine to model real-world systems. The governing wave equation describes the evolution of transverse displacement over time, while modal analysis demonstrates how seemingly complex vibrations may be decomposed into simpler natural modes. Expressing the membrane motion as a combination of its natural modes greatly simplifies the problem while preserving its essential physical behavior. The techniques employed in this investigation extend well beyond membrane dynamics, finding widespread application in structural dynamics, aerospace engineering and numerous other areas of computational physics. The djembe therefore serves as a practical example through which several important computational physics techniques can be applied and studied.

Previous Research

The theoretical foundations of vibrating membranes were established through the classical work of Rayleigh (1877) and later summarized within the broader theory of structural vibrations by Timoshenko (1955), who analysed the vibration of circular membranes under ideal boundary conditions. These works, which used separation of variables and Bessel functions, demonstrate that membrane vibrations may be expressed as a superposition of natural modes whose spatial distributions are determined by the membrane geometry and whose temporal behaviour is governed by harmonic oscillation.
​
Although these analytical solutions remain fundamental to the theoretical understanding of membrane dynamics, they are generally restricted to idealized geometries and boundary conditions. To overcome these limitations, previous computational studies have extended classical membrane theory using numerical methods to simulate the behaviour of percussion instruments under more realistic conditions. Lima-Rodriguez et al. (2019), for example, used the Finite Element Method to compare numerically computed natural frequencies against analytical membrane theory across a range of tension regimes. Other studies have incorporated damping and coupling with surrounding structures directly within a finite difference time-stepping scheme: Bilbao (2010) simulates the coupling between two membranes and an intervening set of snare wires, while Torin (2015) develops a broader finite difference framework for three-dimensional percussion instrument synthesis, including membrane and air-cavity coupling. In many finite element models, by contrast, the membrane displacement is instead represented using modal decomposition, allowing the governing partial differential equation to be transformed into a finite set of coupled ordinary differential equations describing the temporal evolution of each vibration mode; Avanzini & Marogna (2010) apply this modal approach to the sound synthesis of membrane percussion instruments.
​
Research within musical acoustics has also investigated impact mechanics, contact modelling and acoustic radiation in percussion instruments. Many computational models represent strikes using defined initial displacement or velocity distributions, or through modal excitation methods, without modelling the spatial extent of hand contact. Other models, such as that of Avanzini & Marogna (2010), do incorporate finite-duration contact through a dynamically coupled hammer-membrane system based on the Hunt-Crossley contact law, in which the contact duration emerges implicitly from the impact dynamics rather than being explicitly prescribed. However, this and similar models typically treat contact as occurring at a single point, rather than distributed across a realistic contact area, and are not tailored to differentiate the range of playing techniques characteristic of hand percussion. In contrast, this investigation combines a spatially distributed pressure profile — Hertzian-inspired for the finger and slap strikes, annular for the cupped-palm bass strike, and in each case representing the finite contact area of the hand rather than a single point — with a finite-duration temporal loading profile specific to each technique, and calibrates the resulting damping directly against a recorded djembe strike.
​
The studies cited above have significantly advanced the numerical modelling of membrane vibration and percussion acoustics. However, comparatively little attention has been given to modelling finite-duration hand strikes while also investigating how residual membrane vibration influences subsequent strikes during realistic rhythmic performance.

Research Gap and Novelty

The numerical techniques employed throughout this investigation — including finite element analysis and Rayleigh damping—are well established within computational physics and are not novel. Instead, the originality of this work lies in combining these methods to construct a computational model of djembe performance that incorporates realistic strike geometries, finite-duration contact and experimentally calibrated damping and further the modal interactions in a polyrythm.
​
Many simplified membrane models represent impacts through prescribed initial displacement or velocity conditions, effectively assuming that the striking process has already concluded before the simulation begins. In contrast, the present model represents the striking process itself by applying distributed loading profiles over finite contact durations, using contact pressure distributions tailored to the bass, finger and slap techniques. This allows the transfer of force from the hand to the membrane to be modelled explicitly rather than approximated as an instantaneous point impulse.
​
More importantly, the completed model provides a computational framework for investigating how residual membrane vibration influences subsequent strikes during realistic rhythmic performance. Rather than considering impacts in isolation, the model is used to simulate sequences of strikes corresponding to djembe rhythms and quantify how previously excited vibration modes affect the excitation, energy distribution and vibrational response of later strikes. To the author’s knowledge, relatively little previous work has investigated this phenomenon within a physically based computational model of hand percussion.
​
Consequently, the principal contribution of this investigation is not the development of new numerical methods or contact profiles, but the application of established computational techniques to study modal interaction during a realistic djembe performance.

Objectives

The primary objective of this investigation is to develop a computationally efficient yet physically representative model of djembe membrane vibration capable of reproducing the three striking techniques employed during djembe performances.
​
To achieve this objective, the investigation proceeds through several developmental stages. First, the membrane is represented as a circular finite element domain and its natural vibration modes are computed by solving the corresponding eigensystem. These modes form the basis for a modal representation of the membrane displacement, reducing the governing partial differential equation to a finite system of coupled ordinary differential equations describing the temporal evolution of each vibration mode.
​
Subsequently, increasingly realistic physical effects are incorporated into the model. Rayleigh damping models structural energy dissipation, while frequency-dependent air damping accounts for energy losses due to acoustic radiation. The Rayleigh damping coefficients are calibrated using experimentally recorded audio to improve agreement between the simulated and measured vibrational observations. Distributed contact pressure profiles are used to represent realistic contact geometries for bass, finger and slap strikes. Finite contact durations are explicitly represented using time-dependent force functions, allowing the striking process itself to be simulated.
​
Finally, the resulting membrane displacement is reconstructed through modal superposition and visualized using three-dimensional video outputs. The completed model provides both qualitative visualization of membrane deformation and quantitative predictions of modal displacement, vibration amplitude and energy transfer, making it suitable for investigating modal interaction and residual vibration during a real djembe performance.
​
​
The following section presents the numerical implementation of the computational model.

Numerical Implementation

Defining parameters

The numerical implementation begins by defining the physical geometric parameters, material properties and boundary conditions of the membrane before constructing its finite element representation. Together, these parameters define the system and used to establish the governing equations that are solved throughout the remainder of the analysis.​The membrane is modelled as a circular two-dimensional domain, where the thickness is assumed to be negligible compared to its radius. Its motion is governed by the forced two-dimensional wave equation,​ρₛ
∂²u
∂t²
= T∇²u + f(x,y,t),​where u(x,y,t) denotes the transverse membrane displacement, ρₛ is the membrane surface density, T is the membrane tension and f(x,y,t) represents the external force applied during a strike. Rather than solving this partial differential equation analytically, the membrane is discretized using the Finite Element Method (FEM), dividing the continuous domain into a large number of small triangular elements. This discretization transforms the governing equation into a finite system of algebraic equations while preserving the circular geometry and fixed boundary conditions of the membrane. A sufficiently fine mesh is defined to ensure that higher vibration modes are represented accurately without incurring excessive computational load.​The physical parameters governing the membrane, including its radius, tension, surface density and strike location, are defined at this stage. The membrane wave speed is given by​c =
T
ρₛ
​​where T is the membrane tension and ρₛ is the membrane surface density. These parameters define the mechanical behavior of the membrane and establish the physical system upon which all subsequent computations are based. Where possible, parameter values were selected from measurements of the physical instrument or from values reported in the literature. Parameters requiring calibration were subsequently refined using experimentally recorded djembe audio, ensuring that the numerical model remained both physically and acoustically plausible and representative.​The membrane tension and surface density were selected to produce natural frequencies consistent with those measured from the physical instrument while remaining within values reported in the literature. For a uniform circular membrane, these parameters predict a fundamental frequency of approximately 76.8 Hz and a second natural frequency of approximately 122.4 Hz, closely matching the 76.59 Hz and 122.92 Hz frequencies extracted from the recorded djembe strike. The predicted fundamental also falls within the 65–80 Hz range commonly reported for djembe bass tones in literature, supporting the physical plausibility of the selected parameter values.
Load the Finite Element Method package and set the working directory.
Needs["NDSolve`FEM`"]
Define parameters.
In[]:=
radius=0.15;​​mesh=ToElementMesh[Disk[{0,0},radius],MaxCellMeasure->0.00002];​​x0=0.05;​​y0=-0.05;​​tension=725;​​surfaceDensity=0.8;​​realC=Sqrt[tension/surfaceDensity];​​strikeMagnitude=-1.0;

Determining the natural modes

Having established the computational domain of the model, the membrane’s natural modes of vibration are determined by solving the corresponding eigenvalue problem. Physically, every stretched membrane possesses a discrete set of natural frequencies at which it naturally ‘prefers’ to oscillate. Associated with each natural frequency is a unique mode shape describing the spatial deformation of the membrane during that vibration.​Following finite element discretization, the governing wave equation is transformed into the generalized eigenvalue problem​Kϕᵢ = λᵢMϕᵢ,​where K is the global stiffness matrix, M is the mass matrix, λᵢ are the eigenvalues and ϕᵢ are the corresponding mode shapes. The natural frequencies are then obtained from​ωᵢ =
λᵢ
.​These mode shapes form an orthogonal basis over the membrane, allowing any displacement to be expressed as the modal expansion​u(x,y,t) = Σᵢ qᵢ(t)ϕᵢ(x,y),​where qᵢ(t) are the time-dependent modal coordinates. This modal representation transforms the original partial differential equation into a finite system of ordinary differential equations, significantly reducing the computational complexity while preserving the physical behavior of the membrane.​The membrane boundary is assumed to remain fixed throughout the simulation — set by the Dirichlet Condition, representing the physical constraint imposed by the drum shell.
NDEigensystem is used to find the eigenvalues and eigenfunctions for the setup.
The computed eigenvalues are first inspected to confirm that they are positive, ordered by increasing magnitude, and either distinct or appropriately repeated, consistent with the expected behavior. Particular attention is given to the appearance of repeated eigenvalues, which arise naturally from the rotational symmetry of a circular membrane and correspond to degenerate mode pairs possessing identical natural frequencies but different spatial orientations.
The corresponding eigenfunctions are then examined. Whereas the eigenvalues determine the membrane’s natural frequencies, the eigenfunctions describe the spatial distribution of displacement associated with each vibration mode. Regions of positive and negative displacement oscillate in opposite directions, while nodal lines present remain stationary throughout the vibration. Together, these mode shapes form an orthogonal basis over the membrane, allowing any subsequent displacement to be represented as a weighted superposition of the individual modes and therefore they are the mathematical foundation for the modal analysis throughout the remainder of the investigation.
To facilitate interpretation of the numerical solution, each mode shape is visualized individually. Visualizing the eigenfunctions provides a qualitative verification of the finite element solution by confirming the expected nodal patterns and symmetry of a circular membrane. As the mode number increases, the membrane deformation becomes progressively more complicated. These visualizations also demonstrate how strike location influences modal excitation, since strikes applied near nodal lines excite particular modes only weakly, whereas strikes at regions of large modal displacement produce much stronger excitation. This behavior forms the physical basis for the later modelling of bass, finger and slap strikes.

Solving the membrane equation

Define ω and the force coefficients.
Following modal decomposition, the original governing partial differential equation is reduced to a finite set of independent second-order ordinary differential equations of the form
​
q̈ᵢ(t) + ωᵢ²qᵢ(t) = 0,
​
where qᵢ(t) is the modal coordinate and ωᵢ is the natural angular frequency of the ith mode. Each equation describes the temporal evolution of a single mode and, in the absence of damping and external forcing, undergoes ideal simple harmonic oscillation. Solving these equations analytically establishes the ideal response of the membrane and provides the simplest physically meaningful representation of its vibration.
Solve the resulting ordinary differential equations (ODEs) using the DSolveValue function.
Once the modal solutions have been obtained, the complete membrane displacement is reconstructed through modal superposition,
​
u(x,y,t)=∑ᵢ qᵢ(t)ϕᵢ(x,y),
​
where qᵢ(t) are the time-dependent modal coordinates and ϕᵢ(x,y) are the corresponding mode geometries. At every point on the membrane and every instant in time, the displacement is obtained by summing the contribution of each vibration mode weighted by its modal coordinate. This reconstruction is one of the main advantages of modal analysis, as it accurately reproduces the membrane dynamics while avoiding the need to solve the full partial differential equation directly. An arbitrary point on the membrane is evaluated to verify that the reconstructed displacement behaves as expected before proceeding to visualization.
The Sum function is used to obtain the membrane displacement.

Visualizing resultant oscillations

The reconstructed displacement is then visualized using a three-dimensional plot of the membrane, exported as a .mp4 file. Successive time steps illustrate the propagation of the resulting transverse waves across the membrane while preserving the fixed boundary condition imposed by the drum shell. These visualizations provide a qualitative verification of the numerical solution by confirming the expected oscillatory behavior, nodal patterns and modal superposition predicted by the computational model.
Note: The initial-condition formulation specifies only the relative modal excitation through the eigenfunction values at the strike location. Since no physical impact force or momentum transfer is modelled, the absolute amplitude (here greater than the defined radius of the membrane — which is physically inconsistent) of the resulting vibration is arbitrary. Consequently, this preliminary model is used only to illustrate modal decomposition and undamped free vibration, while physically meaningful vibration amplitudes are obtained later using the force-based contact model. The same applies for the Gaussian velocity model described below.
Note: All visualizations in this notebook are pre-generated by exporting animation frames as an MP4 video, uploading the video to YouTube, and embedding it within this notebook.

Damping

While the undamped model captures the fundamental oscillatory behaviour of the membrane, it assumes that no energy is lost during vibration and therefore predicts perpetual oscillation. Real drum membranes continually dissipate energy through internal material friction and interaction with the surrounding air, causing vibration amplitudes to decay over time. To account for these effects, damping mechanisms are incorporated into the modal equations.
​
Structural damping is approximated using Rayleigh damping, in which the damping matrix is expressed as
​
C = αM + βK,
​
where M is the mass matrix, K is the stiffness matrix, and α and β are the Rayleigh damping coefficients. Rather than prescribing damping independently for every vibration mode, proportional damping assumes that the modal damping characteristics arise from a combination of the system’s mass and stiffness properties. Consequently, only two coefficients are required to determine the damping behavior of the entire system. These coefficients are obtained by specifying target damping ratios for two selected vibration modes, after which the damping characteristics of all remaining modes are determined automatically by calculating α and β, without the need to explicitly calibrate M and K. This approach captures realistic modal decay while remaining computationally efficient.
​
In addition to structural damping, a frequency-dependent air damping model may be incorporated to represent energy losses associated with acoustic radiation and viscous interaction with the surrounding air. Since higher-frequency modes generally radiate sound more efficiently, they experience greater energy loss than lower-frequency modes. Although empirical, this model improves the physical realism of the simulation while maintaining a simple and computationally efficient formulation.
Illustrative model combining Rayleigh damping and frequency-dependent air damping, with ζ1 and ζ2 being used to solve for α and β.
The Rayleigh damping coefficients and the frequency-dependent air damping model presented above are included as illustrative examples of commonly used damping formulations. Their parameter values are not intended to represent the final djembe model and would normally require calibration using experimental data. In the absence of measured damping ratios across multiple vibration modes, selecting these parameters remains uncertain and may significantly influence the simulated decay behavior.
​
Consequently, the this investigation develops a model employing only Rayleigh damping, with its coefficients calibrated using recordings of a real djembe bass strike. This calibration provides physically grounded damping parameters while avoiding the introduction of additional empirical assumptions associated with the illustrative air damping model.

Calibration and zeta and coefficients calculation

Import the bass audio file from Wolfram Cloud.
The recorded audio is first analysed in the frequency domain using a Fast Fourier Transform (FFT). The dominant spectral peaks provide estimates of the membrane’s modal frequencies, which are then used to calibrate the damping model by calculating the ζ values.
Sample rate defined.
A Fast Fourier Transform (FFT) is applied.
BandpassFilter is used to isolate a narrow frequency band around the dominant frequency, allowing the decay in amplitude over time to be measured.
The first several detected peaks were excluded from the exponential decay fit, as they correspond to the initial transient response of the membrane immediately following impact, during which the assumption of simple exponential decay does not yet hold. Similarly, the final peak was excluded because its amplitude approaches the recording’s noise floor, where peak detection becomes increasingly unreliable. Restricting the fit to the region exhibiting clear exponential decay substantially improves the coefficient of determination (R²). The resulting fit, with R² ≈ 0.95, indicates good agreement with the exponential decay model and supports the reliability of the estimated damping ratio, ζ₁. The second dominant frequency, approximately 122 Hz as calculated above, is then used to determine the second damping ratio required for Rayleigh damping calibration.
As with the ζ₁ fit, the earliest peaks are excluded from the ζ₂ fit to avoid the initial transient response, while the final peak is excluded as its amplitude approaches the recording’s noise floor. A slightly larger number of initial peaks is omitted for this frequency band, reflecting the longer transient response before a clear exponential decay is established.
​
In principle, from here, the measured damping ratios could be used to define a linear interpolation model,
​
ζᵢ = ζ₁ + b(ωᵢ − ω₁),
​
to estimate the damping ratios of the remaining vibration modes. However, the measured values of ζ₁ and ζ₂ correspond closely to two of the computational model’s own natural frequencies and therefore provide suitable calibration targets for the standard two-parameter Rayleigh damping model. Consequently, the Rayleigh coefficients are determined directly from these two measured damping ratios rather than introducing an additional interpolation model.
Using the experimentally determined damping ratios, ζ₁ and ζ₂, the Rayleigh damping coefficients, α and β, are now calculated. These coefficients define the proportional damping matrix used throughout the remainder of the simulation, completing the calibration of the damping model.

Gaussian initial velocity

Define σ and the Gaussian initial velocity profile, and calculate the resulting force coefficients.
Having determined the modal excitation coefficients, the resulting system of damped modal equations is solved numerically. At this stage, only Rayleigh damping is included, while frequency-dependent air damping is omitted so that the influence of the calibrated structural damping can be examined independently. Each modal coordinate satisfies
​
q̈ᵢ(t) + (α + βωᵢ²)q̇ᵢ(t) + ωᵢ²qᵢ(t) = 0,
​
subject to the initial conditions
​
qᵢ(0) = 0, q̇ᵢ(0) = Fᵢ,
​
where Fᵢ is the modal excitation coefficient obtained by projecting the Gaussian velocity distribution onto the modal basis. These equations are solved independently for each vibration mode to obtain the damped temporal response of the membrane.
Equations are solved using DSolveValue with the prescribed initial conditions.
Having solved the damped modal equations, the complete membrane displacement is reconstructed through modal superposition using the calibrated, damped modal solutions. The resulting displacement function forms the basis for subsequent visualization.
Displacement field is calculated using the Sum function on qSolutionsDamped.
Next, the reconstructed membrane motion is visualized using a sequence of three-dimensional surface plot frames representing the membrane displacement at successive time intervals. Each frame illustrates the instantaneous deformation of the membrane while preserving the fixed boundary condition imposed by the drum shell. The resulting video file provides a qualitative assessment of the membrane’s oscillatory behaviour, modal interactions and damping characteristics. It must be noted that at this stage, the z-axis scaling is still arbitrary.
The membrane model now incorporates its fundamental physical properties, natural vibration characteristics, and calibrated damping mechanism. The next stage introduces realistic excitation by explicitly modelling the interaction between the hand and the membrane.
​
The initial-condition formulation reproduces the subsequent free vibration of the membrane but assumes that energy is transferred instantaneously. In practice, the player’s hand remains in contact with the membrane for a finite duration, during which force is applied continuously as the membrane deforms.
The excitation is therefore reformulated using a distributed force model that represents both the spatial distribution and temporal evolution of the applied contact force. This approach provides a more realistic representation of percussion mechanics while retaining the computational efficiency of the modal framework.

Force-based formulation

Model a single Hertzian contact profile
Define ε and T0, and model the tension function.
Use NDEigensystem to obtain the eigenvalues and eigenfunctions.
The updated eigenvalues and eigenfunctions are then examined to verify that the modified membrane model remains physically realistic. Introducing a modest radial variation in membrane tension produces small changes in the natural frequencies and corresponding mode shapes while preserving the expected qualitative behavior and symmetry of the circular membrane. These updated modes form the basis for all subsequent simulations and analysis.
Because the radial tension variation is small relative to the mean membrane tension, the ordering of the natural frequencies remains essentially unchanged. Consequently, the second vibration mode continues to remain close to the approximately 122.92 Hz resonance identified from the recorded djembe strike. The natural angular frequencies are now recalculated, after which the Rayleigh damping coefficients are updated to remain consistent with the modified membrane dynamics. These revised parameters are used throughout the remainder of the force-based simulations.
Recalculate ω and hence α and β using the new eigenvalues.

Modelling different types of strikes

Bass strike

The bass strike is modelled as a distributed external force representing the comparatively large area of the palm contacting the membrane near its center. Unlike finger or slap strikes, the resulting pressure is spread over a broad contact region rather than concentrated at a (relatively) smaller contact area. Instead of assuming a uniform pressure distribution, a smooth annular profile is adopted, in which the pressure rises from a small residual value at the center of the palm to a maximum at a defined ring radius before decaying. This reflects the slight cupping of the hand during a bass strike, which brings the outer palm into firmer contact than the hollow center, and departs from the monotonically decreasing Hertzian distribution used for the finger and slap strikes. To ensure that the applied load remains localized, the pressure distribution is truncated beyond a defined contact radius. Although this introduces a small discontinuity at the contact boundary, the pressure at this radius is small and its influence on the resulting modal excitation is negligible, allowing for this approximation.
​
The temporal evolution of the strike is modelled independently of its spatial distribution, or the geometry of how it is applied. The applied force increases smoothly from zero to a maximum before returning continuously to zero over the contact interval. Consequently, the complete external force acting on the membrane is represented as
​
f(x,y,t) = P(x,y)F(t),
​
where P(x,y) denotes the spatial pressure distribution and F(t) is the temporal loading function. This formulation provides a more realistic representation of hand percussion than an idealized instantaneous impulse while avoiding the artificial excitation of high-frequency vibration modes associated with abrupt force application.
We define the contact parameters.
Defining the pressure profile, and temporal and spatial force coefficients for the bass strike. Further obtaining the modal forcing function.
The resulting forced modal equations are then solved numerically. Each modal coordinate satisfies
​
q̈ᵢ(t) + (α + βωᵢ²)q̇ᵢ(t) + ωᵢ²qᵢ(t) = Fᵢ(t),
​
subject to the initial conditions,
​
qᵢ(0) = 0, q̇ᵢ(0) = 0.
​
During the contact interval, the applied force continuously transfers energy into the membrane according to the spatial and temporal loading functions. Once contact ceases, the forcing term becomes zero and each vibration mode evolves as a freely damped oscillator.
​
The complete membrane displacement is then reconstructed through modal superposition using the numerical modal solutions.
Solving the ODEs using NDSolveValue with the prescribed initial conditions.
Evaluating the reconstructed displacement over the finite element mesh yields the full time-dependent membrane motion, providing the basis for obtaining the frames for visualization and further analysis. The reconstructed membrane displacement function is then defined and evaluated at an arbitrary spatial location and time as a simple verification of the numerical implementation before proceeding to visualization.
Reconstructing the membrane oscillations and conducting a check with arbitrary values within the set bounds.

Finger strike

The same framework is now applied to model the finger strike. Since the membrane geometry, material properties, damping model and governing equations remain unchanged, only the excitation parameters are modified. The finger strike is represented by a substantially smaller Hertzian-inspired contact region positioned close to the membrane boundary. Compared with the bass strike, the reduced contact area, shorter contact duration and lower applied force produce a more localized excitation that preferentially transfers energy into higher-order vibration modes.
Defining the excitation parameters for the finger strike.
The finger pressure distribution and corresponding modal forcing coefficients are then computed using the same projection procedure described for the bass strike.
The resulting forced modal equations are solved numerically using the previously established formulation. The complete membrane displacement is subsequently reconstructed through modal superposition, and the displacement function is evaluated at an arbitrary location and time as a simple verification of the numerical implementation.

Slap Strike

The slap strike is modelled using four simultaneous fingertip contacts. Rather than approximating the slap as a single contact area, the pressure distribution is constructed by superposing four Hertzian-inspired profiles. This representation more closely reflects the physical mechanism of a slap strike, in which multiple localized forces collectively excite the membrane. The shorter contact duration and separated contact points excite the membrane differently from the bass and finger strikes.
​
The computational procedure remains unchanged: the pressure distribution is converted into modal forcing coefficients, the forced modal equations are solved numerically, and the membrane displacement is reconstructed through modal superposition.
The combined pressure distribution is then projected onto the modal basis using the same procedure developed for the previous striking techniques, yielding the corresponding modal forcing coefficients.
The forced modal equations are solved numerically using the same formulation as for the bass and finger strikes. The complete membrane displacement is subsequently reconstructed through modal superposition, and the resulting displacement function is evaluated at an arbitrary location and time as a final verification of the numerical implementation.
To investigate modal interaction between repeated strikes of the same type, a second slap strike is introduced on the opposite side of the membrane by mirroring the original slap geometry. This allows the subsequent interaction analysis to isolate the influence of strike timing and membrane dynamics rather than differences between striking techniques.

Generating visualizations

With constructed computational models for the three principal djembe striking techniques used for this investigation, the final stage of the implementation is to generate visualizations of the resulting membrane motion. Separate simulations are performed for the bass, finger and slap strikes using a common output framework, with the displacement scale adjusted for each animation to clearly display the characteristic amplitude of each strike. The bass strike, for example, produces substantially larger membrane displacements than the finger or slap strikes. Because the natural frequencies and calibrated damping coefficients are identical across all three simulations, each vibration mode exhibits the same intrinsic decay behavior regardless of the striking technique. Consequently, a common simulation time window is used, ensuring that differences observed between the animations arise solely from the differing spatial and temporal excitation mechanisms rather than changes to the underlying computational model.
​
The generated animation frames are subsequently exported as video files, providing a clear visual representation of the membrane dynamics throughout each simulation.

Polyrhythms

With established computational models for isolated bass, finger and slap strikes, the investigation is now extended to simulate successive strikes occurring during realistic rhythmic performance. The membrane is allowed to continue vibrating between impacts. Consequently, each new strike occurs on a membrane that may already possess non-zero displacement and velocity, allowing residual vibration and modal interference to influence the subsequent response. This capability represents one of the major advantages of the computational model built and provides the basis for investigating the interaction between successive strikes.

Bass, slap, one-finger

Defining the delayed modal forces for each strike.
Aggregating individual modal forces from all strikes to the total force function.
The combined modal forcing function is substituted directly into the damped modal equations,
​
q̈ᵢ(t)+(α+βωᵢ²)q̇ᵢ(t)+ωᵢ²qᵢ(t) = Fᵢ(t),
​
which are solved over the duration of the sequence of strikes. Any interaction between successive strikes emerges from the numerical solution through the residual vibration remaining from earlier strikes.
Solving the ODEs.
The complete membrane displacement is then reconstructed through modal superposition using the numerical modal solutions, as done earlier, producing the displacement field across the membrane throughout the sequence time.
Reconstructing the membrane displacement field.
The reconstructed displacement field is evaluated over the finite element mesh at successive time steps to generate three-dimensional surface plots. These frames are exported as a video, providing a visualization of the membrane response throughout the rhythm and illustrating how successive strikes interact through the residual vibration of the membrane.

2:3 polyrhythm — one-finger and slap

The per-strike modal forcing functions and hence the total modal forcing functions is computed
The resulting combined forcing function is substituted into the damped modal equations and solved over a total simulation time of three seconds. Although the rhythmic pattern itself occupies only the first two seconds, an additional second of simulation is included after the final strike to allow the residual membrane vibration to decay naturally. This provides a clearer visualization of the damping behavior and any persistent modal interactions generated during the timeframe.
Solve the resulting ODEs and reconstruct the membrane displacement.
The reconstructed membrane displacement is subsequently evaluated over the finite element mesh to generate a sequence of frames, which are exported as a video, providing a visualization of the membrane.

2:3 polyrhythm — slap and slap

The previous 2:3 polyrhythm combined slap and finger strikes, whose differing contact areas, force magnitudes and durations produce inherently different modal excitation patterns. Although representative of a realistic djembe rhythm, this makes it difficult to distinguish effects arising from strike timing from those caused simply by differences between the two striking techniques. To isolate the influence of timing and spatial location, the same computational framework is therefore applied to a second 2:3 polyrhythm consisting of two identical slap strikes on opposite sides of the membrane as defined earlier.
​
One slap follows the two-beat rhythm while the mirrored slap follows the three-beat rhythm. Since both rhythmic layers now employ identical excitation mechanisms, any observed differences in the resulting interaction energy can be attributed to the relative timing and spatial attributes of successive strikes rather than differences in strike type.
The forcing function calculation, ODE solving, membrane reconstruction and export process is repeated using the appropriate variable definitions.

Modal energy interaction analysis

Bass, slap, one-finger

The interaction energy analysis is first applied to the bass–slap–finger rhythmic sequence. To isolate the contribution of each strike, the bass, slap and finger responses are each solved independently. Consequently, each isolated solution differs from the full simulation only through the absence of the remaining strikes, allowing the effects of linear superposition to be examined directly.
Independently solve the ODEs for each kind of strike.
Conduct an order of magnitude super position check by subtracting the sum of individual qSolutions from the aggregate qSolutionsRhythm.
The maximum discrepancy across all modes is on the order of 10^-6, several orders of magnitude smaller than the peak modal energy (~10^-3), consistent with the default numerical tolerances of ‘NDSolve’ rather than any error in the underlying formulation. This confirms that the response is indeed the linear sum of the individual strike responses, and that any nonzero difference in the energy computed below reflects genuine physical interaction rather than numerical artefact.
​
With the linear superposition confirmed, the total vibrational energy — the sum of kinetic and potential energy across all retained modes — is computed at each instant, both for the actual combined rhythm and for the sum of the three isolated strikes. The difference between these two quantities isolates the interaction energy: the portion of the total energy that cannot be explained by the strikes acting independently.
Calculate modal energy and hence total energy by summing up modal energy across modes.
Define actual energy and energy obtained by linearly summing up modal energies. Then, find the difference.
Plot resulting graphs using the Plot function.
For the bass–slap–finger sequence, the interaction energy remained positive throughout, peaking at approximately 2.5×10⁻⁶ at the instant of the slap strike. This corresponds to roughly 0.2% of the peak total energy of 1.2×10⁻³, which is dominated almost entirely by the bass. No notable features appear at the finger strike, since by t = 0.50 s the bass residual has decayed to a small fraction of its initial value and the finger strike is itself the weakest of the three. When decomposed by mode, the fundamental frequency contributed approximately 1.2×10⁻⁷, which is about thirty times the next-largest modal contribution, and was the only mode not to change sign during the sequence. Under the calibrated Rayleigh damping, the modal decay rate scales as α + βωᵢ², so by the time the slap arrives, only the fundamental retains significant amplitude. In this rhythm it is effectively the sole channel through which one strike can influence the next.

2:3 polyrhythm — one-finger and slap

The modal energy analysis is repeated for the mirrored slap–slap polyrhythm using the same procedure developed previously.
Plot the resulting graphs, as before.
The 2:3 finger–slap polyrhythm, in which the two layers coincide at t = 0 and t = 1.0 s, produced interaction of both signs. The interaction energy was small and positive through most of the cycle but fell to two sharp negative minima of approximately –4×10⁻⁸ each, corresponding to around 19% of the cycle’s peak total energy of about 2.1×10⁻⁷. Both minima occur precisely at the coincidence points and at these instants the combined response of the membrane carries measurably less vibrational energy than the two strikes acting independently would deliver. This is because the excitations are partially out of phase across the modes they share, and the cross term subtracts. The modal decomposition departs from that of the previous rhythm. Mode 4 contributed approximately 6×10⁻⁹, exceeding the fundamental’s 2.5×10 ⁹, and modes 2 and 3 were comparable to mode 1 (fundamental). Neither the slap nor the finger strike couples efficiently to the fundamental — both are applied at 0.75 R, far from the central antinode, with contact durations an order of magnitude shorter than the fundamental period — and the interaction is consequently distributed across the higher modes they excite strongly.

2:3 polyrhythm — slap and slap

The modal interaction analysis is repeated for the mirrored slap–slap polyrhythm. The same procedure is followed.
The superposition error here is larger than in the previous two checks, at approximately 4.5 × 10⁻⁹ across all modes against a peak energy of approximately 4.3 × 10⁻⁷. This is expected, as this rhythm contains ten discontinuous strike onsets rather than three. The error remains an order of magnitude below the smallest interaction energies reported, so superposition still holds in this model.
The mirrored slap–slap polyrhythm (2:3), constructed to hold strike type constant and vary only timing and spatial position, inverted this behavior. The coincidences at t = 0 and t = 1.0 s now produced the largest positive interaction energies of the cycle, approximately 8×10⁻⁸ and 7×10⁻⁸ respectively, with the total energy at t = 1.0 s exceeding the linear sum by roughly 20%. Negative interaction appeared instead at t ≈ 0.33 s and t ≈ 1.33 s, where a strike arrives on the opposite side of the membrane a third of a beat out of step. Simultaneity alone therefore does not determine whether the interaction is constructive or destructive . What determines it is whether the modal excitation coefficients of the two strikes agree or oppose in sign, for each mode. The modal decomposition of this rhythm demonstrates the point clearly. Mode 2 contributed approximately 1.5×10⁻¹⁰, two to three orders of magnitude below modes 3 and 4 at 4×10⁻⁸ and 3×10⁻⁸. Mode 2 possesses a single nodal diameter, and the two mirrored slap centers lie symmetrically on either side of it, at positions where the mode shape takes equal and opposite values. Their modal forcing coefficients are correspondingly antisymmetric and the cross term nearly cancels. The degenerate partner mode 3, rotated by 90°, exhibits no such cancellation and dominates instead. This suppression follows from the strike geometry alone and involves no fitted parameter.
These findings address the question posed at the outset: residual membrane vibration measurably influences subsequent strikes. The interaction is dominated by the fundamental where a bass strike is present and by higher modes otherwise, and its sign is set by the relative phase with which the two strikes excite each shared mode.

Results

The completed model was evaluated in three stages: validation of the finite element formulation and damping calibration; simulation of the three principal djembe striking techniques; and analysis of modal energy interaction during successive strikes.

Validation and calibration

The natural frequencies of the uniform-tension membrane were first compared with the analytical solution for an ideal circular membrane. Using the calibrated membrane tension and surface density, the analytical model predicts a fundamental frequency of 76.8 Hz and a second natural frequency of 122.4 Hz, the latter agreeing with the 122.92 Hz resonance identified in the recorded djembe strike to within 0.5%. The predicted fundamental also falls within the 65–80 Hz range reported for djembe bass tones in the percussion acoustics literature, supporting the physical plausibility of the chosen membrane parameters. Calibration of the Rayleigh damping coefficients using the recorded bass strike produced exponential decay fits with R² ≈ 0.95 for both the fundamental and second vibration modes, indicating good agreement with the assumed exponential decay model over the fitted interval.

Strike-specific excitation

The three force-centered excitation models for the three types of strikes produced clearly distinguishable membrane deformation patterns. The bass strike, owing to its broad central contact area, excited predominantly the fundamental vibration mode and produced large-amplitude, slowly evolving deformation. In contrast, the finger and slap strikes have smaller contact regions, shorter contact durations and strike locations nearer the membrane boundary, resulting in the excitation of a greater proportion of higher-order modes and producing finer, more rapidly varying deformation. These results are consistent with the qualitative acoustic distinction between the sustained, low-pitched bass tone and the brighter, more transient finger and slap tones . Before interpreting the interaction results, the linearity of the modal formulation was verified directly. Summing the three individually simulated bass, slap and finger responses reproduced the combined bass–slap–finger simulation to within a maximum discrepancy of approximately 10 ⁻⁶ per mode, compared with peak modal energies of order 10 ⁻³ . This error is approximately three orders of magnitude smaller than the physical quantities of interest and is consistent with the default numerical integration tolerances of Mathematica' s NDSolve rather than any deficiency in the mathematical basis and formulation. Equivalent checks for the two polyrhythms returned errors of order 10 ⁻⁹, again well below the interaction energies. This confirms that the interaction energies discussed below arise from the quadratic definition of vibrational energy rather than numerical error .

Modal interaction between successive strikes

Across all three rhythmic sequences the interaction energy was non-zero and, in every case, at least an order of magnitude larger than the corresponding superposition error. Residual membrane vibration therefore measurably alters the energy transferred by a subsequent strike. Its sign, magnitude and modal distribution, however, vary substantially between rhythms, and the pattern of that variation is the primary result of this study.
​
The specifics for all three strike techniques is discussed in detail above, after their respective energy plots.

Conclusion

This investigation set out to determine whether residual vibration of a djembe membrane measurably influences the energy transferred by a subsequent strike, and if so, under what conditions. To address it, a modal finite element model of a circular membrane was constructed, damped using Rayleigh coefficients calibrated against a recorded bass strike, and excited using spatially distributed, finite-duration force profiles representing the bass, finger and slap techniques. Three rhythmic sequences were then simulated and the interaction energy was computed directly.
​
The answer is affirmative. Interaction energies ranged from 0.2% of peak total energy in the bass-slap-finger sequence to approximately 20% at the coincidence points of a 2:3 polyrhythm (slap-slap). In every case the effect exceeded the numerical superposition error by at least an order of magnitude, and its physical origin was confirmed by verifying that the modal displacements themselves superpose linearly to within integration tolerance. Interaction is a necessary consequence of energy being quadratic in a quantity that adds linearly.
​
The more substantive finding concerns the sign of the interaction. Two identical slap strikes applied simultaneously at mirrored positions reinforce one another, whereas a slap and a finger strike applied simultaneously at different positions partially cancel. Simultaneity is therefore not the governing variable. What governs the interaction is whether the two strikes excite each shared vibration mode in phase or in anti-phase, which is fixed by the geometric relationship between the strike locations and that mode’s nodal characteristics. The suppression of the mode 2 interaction in the mirrored slap–slap rhythm, by two to three orders of magnitude relative to its neighbors, is the clearest demonstration: the two strike centers straddle mode 2’s nodal diameter, their forcing coefficients are equal and opposite, and the cross term vanishes. This emerges straight from the geometry.
​
A third result, less anticipated, is that the mode carrying the interaction is itself rhythm-dependent. Where a bass strike is present, the fundamental dominates: its slow decay under the calibrated damping leaves it the only mode with substantial amplitude when the next strike lands. Where both finger and slap strikes are applied near the rim, with contact durations far shorter than the fundamental period, neither couples efficiently to the fundamental and the interaction shifts to the higher modes. The distinction has a musical corollary. It suggests that what a player controls when placing a strike relative to the preceding strike is not only how much the two interact, but which part of the instrument’s spectrum they interact through.
The investigation therefore closes on a quantitative account of an effect that has previously been described only qualitatively. Closely-spaced djembe strikes do interact; the interaction is measurable, is governed by strike geometry rather than by timing alone, and is channeled through different vibration modes depending on the technique employed. Whether it is also audible is a question this model raises but cannot answer.

Future directions

Several extensions would strengthen and broaden the findings of this investigation.
​
The most pressing is a convergence study in the number of retained modes. Only five were used throughout for computational efficiency, yet the slap–slap polyrhythm places its largest interaction contributions in modes 3 and 4 which is uncomfortably close to the truncation boundary. Until it is shown whether the interaction migrates further upward or converges, no conclusion about its modal distribution can be regarded as secure, despite likely being accurate. A mesh refinement study and a direct comparison of the computed eigenvalues against the analytical Bessel-function solution would similarly strengthen the validation, which currently rests on a single second-mode agreement. Additionally, the modal energies reported here are expressed in the arbitrary units of the modal coordinates, hence ratios are meaningful, but absolute values are not, and mass-normalization of the eigenfunctions would be required before they could be read as physical energies.
​
The damping calibration rests on a single recorded bass strike. Both damping ratios were extracted below 125 Hz, and the Rayleigh damping model then extrapolates across a quadratic in ω to set the decay rates of the higher modes. Those rates are therefore fixed by extrapolation rather than measurement, which matters because the relative persistence of the modes is precisely what governs how strongly one strike overlaps the next. Measuring damping at higher frequencies, and repeating the calibration across multiple strikes and instruments to establish confidence intervals, would place the higher-mode results on the same strong empirical footing as the fundamental. The bass, finger and slap contact parameters — contact radius, force magnitude, contact duration, and for the bass strike additionally the ring radius and center pressure floor — likewise remain estimated rather than measured.
​
The model omits coupling to the enclosed air cavity and radiation damping into the open air, both of which shift modal frequencies and decay rates. Since the investigation concerns modal interaction rather than the radiated sound field, this is a reasonable first-order approximation, but it becomes limiting if quantitative prediction of the sound is required — which is a prospective future direction. A few other simplifications could also be relaxed: the radial tension profile was assigned an illustrative value rather than fitted; the boundary was treated as a rigid clamp using a Dirichlet condition rather than a flexible tensing system; the equations were kept linear, neglecting amplitude-dependent stiffening at large bass-strike amplitudes; and the hand was assumed to leave the membrane instantaneously, whereas in practice a hand left resting on the head adds considerable damping which an effect players exploit deliberately for performance.
​
A further limitation is that the model produces no sound. Its predictions are expressed as modal energies and displacement fields, which can be plotted but not heard, and an investigation motivated by the audible differences between djembe strikes cannot be intuitively understood in those terms alone. Rendering the simulated membrane velocity as a synthesized audio file would allow the predicted interaction to be listened to directly and compared against recordings, connecting the numerical results back to the musical phenomenon that motivated them.
​
Finally, the analysis considers only a three-strike sequence and two 2:3 polyrhythms. Extending the framework to longer, more rhythmically realistic patterns while incorporating the refinements above would move it closer to a quantitatively validated acoustic model of djembe performance.

References

• Avanzini, F., & Marogna, R. (2010). A modular physically based approach to the sound synthesis of membrane percussion instruments. IEEE Transactions on Audio, Speech, and Language Processing, 18(4), 891–902. https://doi.org/10.1109/TASL.2009.2036903
• Bilbao, S. (2009). Numerical Sound Synthesis: Finite Difference Schemes and Simulation in Musical Acoustics. Wiley.
• Bilbao, S. (2012). Time domain simulation and sound synthesis for the snare drum. Journal of the Acoustical Society of America, 131(1), 914–925.
• Lima-Rodriguez, A., Gonzalez-Herrera, A., & Garcia-Manrique, J. (2019). Study of the dynamic behaviour of circular membranes with low tension. Applied Sciences, 9(21), 4716. https://doi.org/10.3390/app9214716
• Marogna, R., Avanzini, F., & Bank, B. (2010). Energy based synthesis of tension modulation in membranes. Proceedings of the 13th International Conference on Digital Audio Effects (DAFx-10), 102–108.
• Rayleigh, J. W. S. (1877). The Theory of Sound (Vol. 1). Macmillan and Co.
• Rhaouti, L., Chaigne, A., & Joly, P. (1999). Time-domain modeling and numerical simulation of a kettledrum. Journal of the Acoustical Society of America, 105(6), 3545–3562.
• Timoshenko, S., & Young, D. H. (1955). Vibration Problems in Engineering (3rd ed.). D. Van Nostrand.
• Torin, A. (2016). Percussion instrument modelling in 3D: Sound synthesis through time domain numerical simulation.

Acknowledgements

I would like to express my sincere gratitude to my mentor, Piero Sanchez, for his unwavering guidance and support throughout this project. I am also grateful to my mentor group peers, Isabelle Archer and Ellie Lee, for their encouragement, camaraderie, and the many enjoyable moments we shared.
​
I extend my appreciation to the program directors — Rory Foulger, Eryn Gillam, Megan Davis, and Cyrus Taylor — for their continued encouragement throughout the program. I also acknowledge the dedication and assistance of the teaching assistants, whose efforts were invaluable.
​
I would also like to express my sincere thanks to Stephen Wolfram for providing me with the opportunity to work on this project topic and for fostering such a constructive learning environment.
​
I would like to thank my advocate, Sebastian Rodriguez, for his initial guidance and for helping me secure this project topic.

CITE THIS NOTEBOOK

A finite element study of modal interaction between successive strikes on a djembe membrane​
by Tanish Kulkarni​
Wolfram Community, STAFF PICKS, July 9, 2026
​https://community.wolfram.com/groups/-/m/t/3753028