Judy Najnudel
email: judy.najnudel@ircam.fr

David Roze
email: david.roze@ircam.fr

Simulation of the Ondes Martenot ribbon-controlled oscillator using energy-balanced modelling of nonlinear time-varying electronic components

The Onde Martenot is a classic electronic musical instrument. This paper focuses on the power-balanced simulation of its ribbon-controlled oscillator, composed of linear, nonlinear as well as time-varying components. To this end, the proposed approach consists in formulating the circuit as a Port-Hamiltonian System, for which power-balanced numerical methods are available. A specificity of the Martenot oscillator is to involve parallel capacitors, one of them having a capacitance which non-linearly depends on the time-varying ribbon position state. In the case of linear time-invariant (LTI) capacitors in parallel, an equivalent component can be deduced using the classic impedance approach. Such a reformulation into a single equivalent component is required to derive a state-space Port-Hamiltonian representation of a circuit. One technical result of this paper is to propose a method to determine such an equivalent component in the non LTI case. This method is applied to the present Martenot oscillator. Then, power-balanced numerical experiments are presented for several configurations: fixed ribbon position, realistic and over-speed movements. These results are examined and interpreted from both the electronic and mechanical points of view.

Introduction

As the audio industry is moving towards the digital era, the question of the preservation of analog machines and instruments is paramount. This question is especially relevant for the Onde Martenot, one of the first electronic musical instruments [START_REF] Leipp | Les Ondes Martenot[END_REF] invented in 1928, for it is no longer produced and some of its components are now obsolete. A satisfying solution consists in modelling its circuit in order to build a virtual instrument, so that the community of composers, musicians and musicologists may at least have access to facsimile. To model electronic circuits for audio applications, the state-space form known as Port-Hamiltonian Systems (PHS) has proven to be a powerful approach as it guarantees the power balance of the considered system, therefore preserves the passivity of simulations [START_REF] Falaize | Modélisation, simulation, génération de code et correction de systèmes multi-physiques audios: Approche par réseau de composants et formulation Hamiltonienne à Ports[END_REF] even when its components are not linear. It is multi-physical (a system can be electrical, mechanical, thermal or a mix as well) and modular (a system made of several connected PHS is still a PHS). Yet for some circuit configurations, a direct statespace form cannot be derived -the circuit is said not to be realizable -and an equivalent circuit must be computed in order to perform simulations. This is the case with parallel capacitors which must be replaced by a single equivalent capacitor. However, when the components involved in the circuit are not LTI, as some are in the controllable oscillator of the Onde Martenot, the classic impedance approach is no longer suitable and equivalent components must be computed through a specifically designed method. This paper is structured as follows: the Martenot controllable oscillator circuit is presented in section 2, with a par-ticular attention drawn to the realization problem it poses. In section 3, the PHS formalism is briefly described. In section 4, a method to compute equivalent components in this formalism is developed. A modelling of the complete oscillator is then derived, and several configurations are simulated in section 5. Finally the simulation results are discussed in section 6.

1 Ondes Martenot ribbon-controlled oscillator and problem statement

Circuit overview

The Onde Martenot, invented by Maurice Martenot in 1928, is one of the first electronic musical instruments and is based on heterodyne processing. Heterodyning is a technique used to shift high frequency signals into the audio domain. In the Onde Martenot specifically, each one of two oscillators generate a high frequency quasi-sinusoidal voltage (around 80 kHz); one is fixed, and the player controls the second frequency using a sliding ribbon. The sum of these two voltages is an amplitude-modulated signal. Its envelope is detected using a triode vacuum tube, producing an audible sound, for which the frequency is the difference between the two oscillators frequencies. The triode vacuum tube in the detector is a nonlinear component which adds harmonics to the signal. This enriched signal is then routed towards special kinds of loudspeakers (called diffuseurs) selected by the musician, adding another layer of coloration to the sound. The oscillators are made of an LC circuit cou- pled to a triode vacuum tube (for amplification), through a transformer. In the controllable oscillator, one of the capacitors of the LC circuit is variable and controlled by the ribbon. Fig. 1 shows that the total capacitor is in fact made of several capacitors connected in parallel, some of them LTI, but one of them time-varying.

Problem statement

Two capacitors C A and C B connected in parallel are equivalent to a single capacitor C C (fig. 2). In the case of a Linear Time-invariant system, the notion of impedance allows to determine the equivalent capacitor. Indeed, denoting the capacitors impedances

Z A = 1 jC A ω and Z B = 1 jC B ω respec- tively, Kirchoff's laws i C = i A + i B and v C = v A = v B yield the relation 1 Z C = 1 Z A + 1 Z B (1)
This relation characterizes entirely the equivalent component C C and gives the value of its capacitance:

jC C ω = jC A ω + jC B ω ⇒ C C = C A + C B (2)
However, this classic impedance approach is no longer suitable for non LTI systems: if we were to naively define impedance by the ratio v/i (transfer function), that of nonlinear capacitors would still depend on the charge q, itself time-dependent. By definition, time-varying capacitors would also yield a time-dependent transfer function. As the ribbon-controlled capacitance depends on the ribbon position, which itself depends on time, an adapted method to characterize the equivalent capacitor is needed. The PHS formalism allows to represent an energy-storing component by its energy function instead of its impedance. We thus may rely on this notion in a non LTI case, as it is more general.

Port-Hamiltonian Systems: formalism and examples

This section recalls basics on port-Hamiltonian systems (PHS). For detailed presentation, readers can refer to [START_REF] Duindam | Modeling and control of complex physical systems: the port-Hamiltonian approach[END_REF] and [START_REF] Van Der Schaft | Port-Hamiltonian systems theory: An introductory overview[END_REF].

Formalism

Here we rely on a differential-algebraic form adapted to multi-physical systems [START_REF] Hélie | Systèmes Hamiltoniens à Ports avec approche par composants pour la simulation à passivité garantie de problèmes conservatifs et dissipatifs[END_REF] [START_REF] Falaize | Passive guaranteed simulation of analog audio circuits: A port-Hamiltonian approach[END_REF], which allows to represent a dynamical system as a network of storage components with their state variable x and total energy of the state H(x), dissipative components with their variable w and constitutive law z(w), and connection ports as control inputs u and their associated outputs y such as u y is the external power brought to the system. The variables are generally time-dependent and can be vectors. If such a system is realizable [START_REF] Brockett | Finite dimensional linear systems[END_REF] [8], the flows and efforts exchanges between the system components are coupled through a skew-symmetric

The skew-symmetry of S guarantees that the system remains passive, i.e there is no spontaneous creation of energy. Indeed, from Eq. ( 3), the scalar product of the efforts and flows yields

E F = E SE = (E SE ) = -E SE = -E F = 0, (4)
meaning that the following power-balance is satisfied

dE dt ∇H(x) dx dt = P ext u y -P diss z(w) w ≥ 0 , (5)
where E = H(x) is the energy, P ext is the (incoming) external power and P diss ≥ 0 the dissipated power. Appendix A.2 and [2] describe a numerical scheme preserving those properties in discrete time.

Capacitors constitutive laws

For linear time-invariant capacitors, the charge q and the voltage v are mapped according to a constitutive law q = Cv, which depends on a unique characteristic constant (capacitance C in Farad). The electric power P = i v received by such a component makes its stored energy E vary as dE dt = P. With current i = dq dt , voltage v = q C and assuming a zero energy for a discharged component, a time integration yields E = H(q) with H(q) = q 2 2C . This energy is sometimes expressed independently of value C as E = qv 2 . For nonlinear capacitors, the last expression in no longer true. But a description based on an energy function q → H(q) is still applicable. The constitutive law is described by the voltage function H (derivative of H), namely,

v = H (q), (6)
with power balance v i = H (q) dq dt = dH(q) dt = dE dt .

Remark 1 (Constitutive laws based on q or v). In practice, constitutive laws are usually formulated (and measured) with respect to the voltage (effort) rather than the charge (state). Formally (if possible), such a description corresponds to invert v = H (q) (= q/C for linear capacitors) into q = F(v) (= Cv for linear capacitors) with F = H -1 and formulates the energy as

E = H H -1 (v) = Cv 2 2 .
Note that differentiating this formula yields a (correct) power balance but more difficult interpretations. This remark applies to varactors proposed in [START_REF] Sarti | Toward nonlinear wave digital filters[END_REF], with model

H -1 = F : v → Cv/ 1 + v/v 2 .
Another example is tanh type which is of the form q → H (q) = v 1 tanh( q q 0 ). Fig. 3 shows the voltage functions and the corresponding energy functions of these different capacitors types (for

q 0 = 7nC, v 1 = 80mV, C = 50nF and v 2 = 0.2V).

Equivalent component description of non-LTI parallel capacitors

Method 3.2 Problem statement and hypotheses

Consider two capacitors connected in parallel (fig. 2). These components are flow-controlled ( dx dt = q = i, current). As they are connected in parallel, their dual efforts (voltages) are equal and there is no skew-symmetric matrix S such that

  i A i B -V   = S.   v A v B I   ,
and formulation (3) cannot be retrieved. Replacing those parallel capacitors by a single equivalent capacitor (with common voltage v C := v A = v B and total current i C = i A + i B ) allows to restore such a formulation. For linear time-invariant (LTI) laws, the use of transfer functions of impedance-type makes this operation straightforward. For non LTI laws, no characterization can be based on usual transfer functions. In this case, energy-storing components are characterized by their energy function to be used in the PHS formalism.

The purpose of this section is to derive the energy function H C of the equivalent component from the energy functions H k (k = A, B) of isolated components, including for non LTI laws, under the following hypotheses:

(i) the energy function H k is C 1 positive-definite (H k (0) = 0 and H k (x) > 0 for x = 0) (ii) the voltage function H k (derivative of H k ) is strictly increasing and definite (H k (0) = 0)
According to Eq. ( 6), this means that the voltage v k = H k (q k ), continuously and strictly increases with the charge q k , and is zero for a zero charge. In particular, the constitutive law Eq. ( 6) makes v k and q k in one-to-one relation, allowing therefore its invertibility.

Method

To express the total energy H C as a function of the total charge q C = q A + q B under the constraint that capacitors in parallel share the same voltage v A = v B = v C , the method is decomposed into three steps.

Step 1: express the total charge q C = q A + q B as functions of the common voltage

v A = v B = v C .
From Eq. ( 6), the charge of isolated components is

q k = H k -1 (v k ), for k = A, B, (7)
so that the total charge depends on the common voltage as

q C = H A -1 + H B -1 (v C ). (8)
This function continuously and strictly increases and is zero at zero.

Step 2: express this common voltage v C as a function of the total charge q C .

v C = H A -1 + H B -1 -1 (q C ). ( 9
)
Step 3: express the total energy as a function of q C . The energy values H k (q k ) of elementary components k = A, B can be reformulated as functions of the total charge, using the composed functions q C (9)

-→ v C (7)
-→ q k . Their sum yields the total energy function, that is,

H C (q C ) = H A • H A -1 + H B • H B -1 • H A -1 + H B -1 ] -1 (q C ). ( 10
)
These steps are detailed on examples in appendix A.3.

Remark 2 (Time-varying case). For capacitors that depend on other additional state variables (e.g. the time-varying space variable in Eq. 15, section 4.1.1), steps 1 to 3 are unchanged (these additional variables are considered as parameters in this method). This is applied to the ribbon controlled-oscillator in section 4.

Generalizations

This method can also be extended to K non-LTI capacitors connected in parallel, leading to

H tot (q tot ) = K ∑ k=1 H k • H k -1 • K ∑ k=1 H k -1 -1 (q tot ) (11)
In this case, the charge of each component k is

q k = H k -1 (v C ) = H k -1 • K ∑ k=1 H k -1 -1 (q tot ) (12)
This method is adapted to other types of storage components that poses similar realization problems. For instance, it is suitable for coils in series, for which the state is the magnetic flux φ , i = H (φ ) provides the current law and the voltage is v = φ . Moreover, when the constitutive laws of the components are not well known, this method can still be used with laws interpolated from measurements. An implementation of the method using piecewise linear functions is available in the PyPHS library [START_REF] Falaize | PyPHS: Passive modeling and simulation in python[END_REF].

4 Simulation of the ribbon-controlled oscillator

Ribbon-controlled oscillator modelling

Constitutive laws of components R 1 , L, C 1 , C 2 , C 13 and C 14 are supposed to be linear and described in table 4, in appendix A.1. The transformer is also supposed to be linear, and of ratio M. The remaining components models are described in the next section.

Equivalent variable capacitor

The ribbon slides between the faces of a comb-shaped capacitor. As it slides, its conductive part hides and activates complementary parts of the teeth that compose the capacitor, according to the ribbon position (fig. 4). Section 4.2 validates that the LTI capacitors C 13 and C 14 and the variable capacitor C 15 connected in parallel are equivalent to a single variable capacitor. In order to evaluate how the equivalent capacitor behaves with respect to its charge q and the ribbon position x, the position is mapped with the heard frequency. The ribbon follows a dummy piano keyboard on which the width

x 0 = 11.10 -3 m ( 13
)
is a semitone [START_REF] Courrier | Analyse de fonctionnement Onde 169[END_REF].

Fig. 4. Variable capacitor of the Onde Martenot with its ribbon control. h is the ribbon height, l is the width of a tooth, x is the ribbon position and d the distance between the ribbon and the capacitor.

PAPERS

Simulation of the Ondes Martenot ribbon-controlled oscillator using energy-balanced modelling of nonlinear time-varying electronic components

Denoting n the semitone number where the reference is A 1 (110 Hz, n = 0), we roughly have x = nx 0 (except when the ribbon slides between a B and a C or between an E and a F where the displacement is greater). The Onde Martenot is tuned on equal temperament, therefore the heard frequency from A 1 is f m = A 1 2 x 12x 0 . Denoting F the carrier frequency, the actual oscillator frequency is f = Ff m and the corresponding capacitance of the LC circuit for a static configuration is (ignoring the dissipative effects of the triode vacuum tube for simplicity):

C(x) = 1 4π 2 (F -A 1 2 x 12x 0 ) 2 L . ( 14
)
Fig. 5 shows C in function of x with L = 2mH and F = 10kHz. Based on the linear electric behaviour v = q C(x) observed for any static position x, the total energy of the electromechanical component has the form H(q, x) = H(q = 0, x) + q 0 ξ C(x) dξ . Moreover, the discharged component (then purely mechanical) applies no force F on the ribbon whatever the position x. This means that F = ∂ H ∂ x is zero at any state (q = 0, x), leading to H(q = 0, x) = H(q = 0, x = 0) = 0: the discharged component stores no pure mechanical energy. Finally, the internal energy of the electromechanichal component is given by (see also fig. 6)

H cap (q, x) = q 2 2C(x) . ( 15
)
Remark 3 (Energy time-variation). The total energy E = H(q, x) varies as dE dt = P e + P m where r P e = ∂ H ∂ q (q, x) dq dt accounts for the incoming electrical power due to current dq dt and voltage v = ∂ H ∂ q (q, x) = q C(x) , r P m = ∂ H ∂ x (q, x) dx dt accounts for the incoming mechanical power due to velocity dx dt and a spring reaction force F = ∂ H ∂ x (q, x) = -C (x)q 2 2C(x) 2 , induced by the capacitance variation.

Triode

The triode vacuum tube is modelled with an enhanced Norman Koren model [START_REF] Cohen | Modélisation, analyse et identification de circuits non linéaires: Application aux amplificateurs guitare à lampes pour la simulation en temps réel[END_REF]. This gives the anode current i pc and grid current i gc in function of the voltages v pc and v gc :

i pc = 2E E x 1 /K g if E 1 ≥ 0 0 else (16) i gc = 0 if v gc < V a v gc -V a R gk else ( 17
)
with

E 1 = v pc K p ln(1 + exp(K p ( 1
µ + v gc + V ct K vb + v 2 pc )))
The parameters set θ = (µ, E x , K g , K p , K vb ,V ct ,V a , R gk ) is retrieved from the datasheet [START_REF] Tung-Sol | 6F5 typical operating conditions and characteristics[END_REF] through a least squares minimization. This allows the modelling of the triode as a dissipative component in the PHS formalism, with w = v pc v gc and z θ (w) = i pc i gc .

Transformer and feedback loop

The oscillation starts with a voltage noise V start at the triode grid, which is amplified and filtered by the LC circuit at the triode plate before being re-injected in the grid through a transformer of ratio M. Under open circuit (negligible influence of the following stages) and small-signal (possible linearization around an operating point) hypotheses, the oscillator can be represented as a system of input V start and output V out (fig. 8, with µ the triode amplification factor and ra the anode resistor). The plate load is constituted by a LC parallel circuit, with impedance (s is the Laplace variable)

Z = sL 1 + s 2 LC (18)
Moreover, as V bias is constant, ∆V out = -Z/(Z + ra) ∆v pc where operator ∆ denotes variations around the operating point. Considering the triode amplification factor definition, this yields ∆V out = Z/(Z + ra) µ(M∆V out + V start ), so that

∆V out V start = Zµ Z(1 -µM) + ra . ( 19
)
The poles of transfer function ( 19) are found to be the roots (in s) of the characteristic equation

s 2 + s 1 -µM raC + ω 2 0 = 0 ( 20
)
with ω 2 0 = 1 LC . The condition for the system to start oscillating is that a complex pole has a positive real part. This leads to 1

µ ≤ M < 1 µ + 2ω 0 C gm ( 21
)
where gm = µ/ra is the triode transconductance. The closer to 1/µ is chosen M, the more stable is the oscillation amplitude.

Fig. 8. Schematic of the simplified oscillator.

Interconnection

Replacing parallel capacitors C 1 //C 2 with equivalent capacitor C 3 and C 13 //C 14 //C 15 with equivalent capacitor C 8 , the oscillator is realizable and can be represented as the following PHS:

            v L i C8 i C3 v R1 v pc v gc -I bias -I start             =             0 1 0 0 0 0 0 0 -1 0 0 0 1 -M 0 0 0 0 0 -1 1 1 0 0 0 0 1 0 0 0 0 0 0 -1 -1 0 0 0 1 0 0 M -1 0 0 0 0 1 0 0 0 0 -1 0 0 0 0 0 0 0 0 -1 0 0             .             i L v C8 v C3 i R1 i pc i gc V bias V start            

Numerical experiments

We now simulate the complete variable oscillator, sliding the ribbon from A1 to A6. The simulation is performed ac-cording to the power-balanced numerical scheme presented in [START_REF] Falaize | Modélisation, simulation, génération de code et correction de systèmes multi-physiques audios: Approche par réseau de composants et formulation Hamiltonienne à Ports[END_REF] (see also the PyPHS library [START_REF] Falaize | PyPHS: Passive modeling and simulation in python[END_REF]). In order to observe the frequency changes, we choose a carrier frequency F = 10kHz (instead of the actual 80 kHz) and perform the sweep in 2 ms (vel 1 ) or in 0.5s (vel 2 ). The simulation parameters are presented in table 1. Due to important nonlinearities of some components, a large sample-rate is chosen to avoid any aliasing.

V ct V a R gk 0.5 0.33 1300 R R R 1 1 1 L L L C C C 1 1 1 //C C C 2 2 2 C C C 13 //C C C 14 V V V bias M M M 7.5kΩ 2mH 0.22µF 440pF 90V 1 µ
Table 2. Components parameters values. [START_REF] Falaize | PyPHS: Passive modeling and simulation in python[END_REF] show the observed flows and efforts of the oscillator as the ribbon slides. Figs. 11 and12 show the states of the different storage components, reflecting the frequency changes. Figs. 14 and 15 show the power balance of the complete system during the simulation. Fig. 13 shows the spectrogram of the output voltage, suggesting that the harmonic distortion is sufficiently negligible. It is also worth noting that the PHS formalism gives access to other physical parameters of the ribbon, which would otherwise prove difficult to measure. Indeed, the quantity

∇ x H cap (x, q) = -q 2 C (x) 2C 2 (x) (22)
is the force produced by the ribbon displacement (see Remark 3) ; figs. 16 and 17 show that the values taken by this force during the sweep are negligible no matter the sweep velocity.

q q qC C C 3 3 3 r r re e e f f f q q qC C C 3. Reference values.

Interpretation and discussion

Due to the nature of the Ondes Martenot instrument (rare, fragile and expensive), setting up extensive measurements to evaluate the accuracy of the oscillator model is a complex operation which is still ongoing. However, a preliminary observation is that with the chosen parameters, the oscillation is quasi-sinusoidal (less that 0.1 % harmonic distortion for the second harmonic), which corresponds to observations made in [START_REF] Courrier | Analyse de fonctionnement Onde 15[END_REF] : "on the whole ribbon range, the sinusoidal quality of the signal produced by the oscillator is excellent". Moreover, the behaviours of the components are not affected by the ribbon displacement speed and there is no latency between the ribbon displacement and the frequency changes. This suggests that for future simulations of the complete circuit, the oscillator could be modelled with a frequency-controlled sinusoidal voltage generator, which would save computation time. It can be noted though that the output voltage of the oscillator slightly decreases as 3 shows the reference values. 3 shows the reference values. the heard frequency increases, which is not the case in real Ondes Martenot. Real instruments have a higher frequency carrier, therefore frequency modulations are relatively less important. The oscillator capacitance varies less and consequently the triode amplification, which depends on its load, is stable during the whole sweep. A second observation is that consistently with the PHS formalism, the power balance is preserved despite the nonlinearities of some components. A third observation is that the mechanical force produced by the ribbon displacement is extremely low (less than 10nN), thus presumably not detectable by the musician. This is likely to be what Maurice Martenot intended as the interaction between the player and the instrument should be as smooth as possible.

Conclusion and perspectives

In this paper, a refined simulation of the ribbon-controlled oscillator of the Ondes Martenot has been proposed, allowing a numerical investigation of this circuit. It is based on an energy-balanced modelling adapted to LTI (capacitors, resistor, inductor) and non LTI (vacuum tube and the multiphysical time-varying capacitor mechanically-driven by a ribbon) components. One contribution of this paper is the design of a method to compute equivalent energy functions of groups of components when required, to derive a PHS statespace realization of a global circuit. This method reveals that LTI and variable capacitors in parallel build an equiva- lent bi-variate capacitor, depending on an electric state (total charge) and a space state (ribbon position). As a second contribution, the ideal energy function for the ribbon-controlled equivalent capacitor is derived in correspondence with the target keyboard designed by Martenot. Numerical experiments on the nonlinear time-varying circuit lead to expected observations: (1) the combination of the triode amplification and the LC-resonator produces a quasi-sinusoidal oscillation with a stable amplitude for a static configuration; (2) the mechanical force produced by the variable capacitor due to the ribbon displacement is undetectable by the musician (less than 10nN) for over-speed movement ( 300 m/s); [START_REF] Duindam | Modeling and control of complex physical systems: the port-Hamiltonian approach[END_REF] the latency between the instantaneous frequency and the ribbon position is also undetectable. This corroborates that the Martenot's ribbon-controlled circuit is close to an ideal oscillator.

Further work aims to model and simulate the complete Martenot instrument, including all its stages and the diffuseurs. Furthermore, for real-time sound synthesis purposes, the examination of anti-aliasing methods [START_REF] Muller | Trajectory anti-aliasing on guaranteed-passive simulation of nonlinear physical systems[END_REF] could be profitable to handle the nonlinear heterodyne process at reasonable sampling rates. pacitor C and the inductor L are storage components whose states are given by the variables q (charge) and φ (magnetic flux) respectively ; the resistor R is a dissipative component described by Ohm's law. The system is current (flow) controlled, its associated output is a voltage (effort). Table 4 recaps the variables and the constitutive laws of the three components. Kirchoff's laws in receptor convention yield

x dx dt H(x) ∇H(x) C q q = i C q 2 2C q C = v C L φ φ = v L φ 2 2L φ L = i L w z(w) R v R v R R = i R
    i C v L v R -V     =     0 -1 -1 1 1 0 0 0 1 0 0 0 -1 0 0 0     .     v C i L i R I    
Currents are flows and voltages are efforts, therefore their products are powers. Eq. ( 5) is thus naturally retrieved. Note that in this simple example, the matrix S is sparse with constant coefficients, but the properties of the PHS formalism hold for nonlinear or coupled systems which yield more complex matrices.

A.2 PHS : numerical scheme for simulations

The PHS formalism guarantees the passivity of the system in continuous time. Introducing discrete gradient in the numerical scheme [START_REF] Itoh | Hamiltonian-conserving discrete canonical equations based on variational difference quotients[END_REF] allows to preserve this passivity property in discrete time, therefore granting the stability of the simulation as well. Here a one-step numerical scheme is used, yielding

x(k + 1) = x(k) + δ x(k) (23)
In the mono-variate storage component case (H(x) = ∑ N n=1 H n (x n ) where N is the storage components number), the discrete gradient [ ∇H(x, δ x)] n is defined by

[ ∇H(x, δ x)] n = H n (x n +δ x n )-H n (x n ) δ x n if δ x n = 0, H n (x n ) else. ( 24
)
The discrete energy variation is retrieved by chain derivation:

δ E(k) t s = ∇H(x, δ x(k)) δ x(k) t s ( 25
)
where t s is the sampling period. The simulation is achieved replacing dx dt with δ x(k) t s and ∇H(x) with ∇H(x, δ x) in (3). This yields a dynamic equation of the form δ x(k) = t s f k (x(k), δ x(k))

(26)

where f k is a function depending on H, z, u and S. Its solving (using Newton-Raphson method for instance) allows to compute x(k + 1) and y(k).

A.3 Equivalent component method : examples

The method steps (sec. (NL) capacitors with homogeneous laws of common degree α > 0 (nonlinear for α = 1, see Fig. 19).

Label Energy function H k

Voltage function

H k LTI H k (q k ) = q 2 k 2C k v k = H k (q k ) = q k C k NL H k (q k ) = E k q k q 0 1+α v k = H k (q k ) = V k q k q 0 α
In the model (NL), [x] α = sign(x) |x| α denotes the signed power function and constant values q 0 (charge), E k (energy) and V k (voltage) are related by

V k = (1 + α)E k /q 0 . ( 27
)
This leads to the following sequence of derivations.

Step 1. Eq. ( 7) yields q k = H k -1 (v k ) with (LTI):

H k -1 (v k ) = C k v k , ( 28
) (NL): H k -1 (v k ) = q 0 v k V k 1/α , (29)
and Eq. ( 8) yields

q C = Q C (v C ) with Q C (v C ) := H A -1 + H B -1 (v C ), (LTI): Q C (v C ) = (C A + C B )v C (30) (NL): Q C (v C ) = q 0 v C V A 1/α + q 0 v C V B 1/α = q 0 v C V C 1/α , (31)
with

V C = V -1/α A + V -1/α B -α
.

Step 2. Eq. ( 9) yields v C = Q -1 C (q C ) with (LTI):

Q -1 C (q C ) = q C C A + C B ( 32
) (NL): Q -1 C (q C ) = V C q C q 0 α . ( 33
)
Step 3. Eq. ( 10) yields (LTI):

H C (q C ) = ( C A q C C A +C B ) 2 2C A + ( C B q C C A +C B ) 2 2C B = q 2 C 2(C A + C B ) (34)
(NL):

H C (q C ) = ∑ k=A,B E k q 0 V C q C q 0 α V k 1/α q 0 1+α = E C q c q 0 1+α ( 35
)
where E C is found to be related to V C as E k to V k in Eq. ( 27). Note that the equivalent laws for (LTI) and (NL) have the same expression. This is not the case in general.

THE AUTHORS

Fig. 1 .

 1 Fig. 1. Schematic of the Onde 169 controllable oscillator (source: Musée de la Musique, Paris).

Fig. 2 .

 2 Fig.2. Equivalence between two parallel capacitors and a single capacitor.

Fig. 3 .

 3 Fig. 3. (a) Voltage functions H (associated with constitutive law Eq. (6)) and (b) energy functions H for different capacitors types.

Fig. 5 .

 5 Fig. 5. Capacitance of the variable oscillator in function of the ribbon position from A 1 to A 6 .

Fig. 6 .Fig. 7 .

 67 Fig. 6. Energy function of the variable capacitor.

8 8 8 r r re e e f f f φ φ φ r r re e e f f f vel 1 4 .

 4 60.10 -6 C 3.13.10 -8 C 4.60.10 -6 Wb vel 2 4.58.10 -6 C 4.76.10 -8 C 4.58.10 -6 Wb Table

Fig. 9 .

 9 Fig. 9. Velocity vel 1 = 384 m/s: simulated flows and efforts of the storage components, for a carrier frequency F = 10kHz, during a sweep from A1 to A6.

Fig. 10 .

 10 Fig. 10. Velocity vel 2 = 1.54 m/s: simulated flows and efforts of the storage components, for a carrier frequency F = 10kHz, during a sweep from A1 to A6 zoomed on the first 2ms.

Fig. 11 .

 11 Fig. 11. Velocity vel 1 = 384 m/s: simulated states of the storage components during a sweep from A1 to A6. Table3shows the reference values.

Fig. 12 .

 12 Fig. 12. Velocity vel 2 = 1.54 m/s: simulated states of the storage components during a sweep from A1 to A6 zoomed on the first 2ms. Table3shows the reference values.

Fig. 13 .

 13 Fig. 13. Spectrogram of the output voltage during a sweep from A1 to A6.

Fig. 14 .Fig. 15 .

 1415 Fig. 14. Velocity vel 1 = 384 m/s: power balance of the system during a sweep from A1 to A6.

Fig. 16 .

 16 Fig. 16. Velocity vel 1 = 384 m/s: mechanical flow and effort during a sweep from A1 to A6.

Fig

  Fig. 18. Parallel RLC

  Fig. 19. (NL) Energy functions H k (top) and voltage functions H k (bottom) for various degrees α (q 0 = 0.9.10 -8 C and V k = 0.5V).

Table 1 .

 1 Table 2 recaps the components values. We denote C 1 //C 2 = C 3 and C 13 //C 14 //C 15 = C 8 . Figs. 9 and Simulation parameters.

		F F Fs s s	F F F	v v ve e el l l 1 1 1	v v ve e el l l 2 2 2
		768kHz 10kHz 385 m/s 1.54 m/s
	6F5	µ	E x	K g	K p	K vb
		98	1.6	2614	905 1.87

Table 4 .

 4 State variables and constitutive laws for a linear parallel RLC circuit.

J. Audio Eng. Sco., Vol. , No. ,

J. Audio Eng. Sco., Vol. , No. ,

Acknowledgments

This work was funded by the Collegium Musicae and is supported by the project ANR-16-CE92-0028, entitled Interconnected Infinite-Dimensional systems for Heterogeneous Media, INFIDHEM, financed by the French National Research Agency (ANR).


Marc Wijnand
email: marc.wijnand@ircam.fr

Brigitte D'andréa-Novel
email: brigitte.dandreanovel@ircam.fr

Benoît Fabre
email: benoit.fabre@upmc.fr

Thomas Hélie
email: thomas.helie@ircam.fr

Lionel Rosier
email: lionel.rosier@mines-paristech.fr

David Roze
email: david.roze@ircam.fr

B D'andréa

Active control of the axisymmetric vibration modes of a tom-tom drum

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Active control of the axisymmetric vibration modes of a tom-tom drum Marc Wijnand, Brigitte d'Andréa-Novel, Benoît Fabre, Thomas Hélie, Lionel Rosier, and David Roze Abstract-This paper deals with an application of active control of percussion instruments. Our setup consists of a tom-tom drum with a circular membrane, a cylindrical cavity and a circular rigid wall on which a loudspeaker is mounted. The current applied to the loudspeaker is controlled in order to modify the frequencies of the drum membrane modes. First, a PDE model of the axisymmetric transverse vibration of the tom-tom membrane is developed. Subsequently, the equation is recast as an infinite-dimensional port-Hamiltonian system. The port-Hamiltonian framework enables us to develop a numerical scheme that preserves the power balance and guarantees a stable simulation. Finally, a control law for the loudspeaker current is designed to modify the frequency of the first axisymmetric vibration mode of the drum membrane, using finite-time and passivity-based methods.

I. INTRODUCTION A. Active control of musical instruments

Active control of musical instruments consists in adding a control loop to an existing acoustic musical instrument that is being played by a musician. It is an example of active vibration control [START_REF] Elliott | Active noise control[END_REF], [START_REF] Fuller | Active control of vibration[END_REF] that enables the musician to enlarge his sound palette while keeping the ergonomics of the original instrument. In addition, it can be used to study the acoustics of the original instrument.

One can distinguish two classes of active control of musical instruments [START_REF] Andréa-Novel | An acoustic model for automatic control of a slide flute[END_REF], [START_REF] Boutin | Méthodes de contrôle actif d'instruments de musique. cas de la lame de xylophone et du violon[END_REF], [START_REF] Benacchio | Contrôle actif modal appliqué aux instruments de musique à cordes[END_REF], [START_REF] Meurisse | Contrôle actif appliqué aux instruments de musique à vent[END_REF], [START_REF] Jossic | Modal active control of Chinese gongs[END_REF]. In acoustical active control, the control acts on a fluid medium. An example is the use of a loudspeaker to create destructive interference in order to cancel sound. In structural acoustical control, the control acts on a solid. Examples include the use of piezoelectric actuators placed on a violin, guitar or percussion instrument.

B. The tom-tom drum

The tom-tom drum is an unpitched directly struck membrane instrument. It consists of a cylindrical body with a top and bottom membrane and is a standard part of a drum kit.

In the field of acoustics of percussion instruments with an air-membrane coupling, instruments such as the tom-tom drum and snare drum [START_REF] Worland | Normal modes of a circular drumhead under nonuniform tension[END_REF] and Indian tabla [START_REF] Sathej | The eigenspectra of Indian musical drums[END_REF], [START_REF] Tiwari | Effects of air loading on the acoustics of an Indian musical drum[END_REF] have been considered. Most research is devoted to the timpani or kettledrum [START_REF] Christian | Effects of air loading on timpani membrane vibrations[END_REF], [START_REF] Rhaouti | Time-domain modeling and numerical simulation of a kettledrum[END_REF].

C. Control type

The drum membrane is modeled as a partial differential equation (PDE) and the loudspeaker as an ordinary differential equation (ODE). The drum membrane and the loudspeaker are coupled by the air in the cavity volume. We can make following remarks about the control.

1) Distributed (modal) control of a PDE: Depending on the positioning of actuators and sensors, two types of control of PDEs can be defined [START_REF] Krstic | Boundary control of PDEs: A course on backstepping designs[END_REF]. In in-domain control, actuators and sensors are located inside the domain. In the case of a lot of distributed actuators and sensors, it is called distributed control. In the alternative case of boundary control, actuators and sensors are located at the boundary of the PDE. In our coupled PDE-ODE system, the PDE is submitted to a distributed pressure that varies as function of a control loop in the ODE.

Furthermore, we opted for a modal control of the PDE. The PDE is projected on eigenfunctions. In the current work, only the first mode is considered.

2) Finite-time control of an ODE: We use finite-time control, a nonlinear control type where the equilibrium point is reached in a finite time [START_REF] Haimo | Finite time controllers[END_REF], which is useful for timeconstrained and robust control [START_REF] Bhat | Finite-time stability of continuous autonomous systems[END_REF], [START_REF] Bernuau | Robust finite-time output feedback stabilisation of the double integrator[END_REF].

II. PHYSICAL MODELLING A. Layout

The tom-tom drum coupled to the cavity and loudspeaker is depicted in Figure 1.

B. Hypotheses

The following hypotheses are made in our model. (H1a) The membrane executes small vertical displacements and its tension is uniform. (H1b) The vibration of the membrane contains only axisymmetric modes 1 . We suppose that the musician excites the membrane at its center. (H1c) The sound radiation in air by the membrane is modeled by adding 0.85aρ air to the membrane's areal density σ [START_REF] Kinsler | Fundamentals of acoustics[END_REF], a being the membrane's radius and ρ air the volumetric density of air. (H2a) The pressure inside the cavity is uniform (lumped (H2b) The air inside the cavity behaves as an ideal gas that is compressed adiabatically. (H3a) The loudspeaker is modeled as a plane piston.

V (t) p(t)
(H3b) The loudspeaker operates above its resonance frequency.

C. Physical model

The model of our system consists of three parts: the membrane (PDE) and the loudspeaker (ODE), that are coupled by the cavity (static relation), as depicted in Figure 1.

1) Membrane: Under hypotheses (H1a)-(H1b), one can state the wave equation with a pressure source as [START_REF] Graff | Wave motion in elastic solids[END_REF] T ∆y m (r, t)

+ p c (t) = σ m ∂ 2 ∂t 2 y m (r, t) + µ ∂ ∂t y m (r, t), (1)
where ∆ = ∂ 2 ∂r 2 + 1 r ∂ ∂r in polar coordinates, y m (r, t) is the transverse displacement of the membrane, T is the tension [N/m], σ m = σ + 0.85aρ air is the equivalent areal density [kg/m 2 ] according to (H1c), µ is a friction coefficient kg/m 2 s and p c is the net pressure [Pa] applied to the membrane. The wave velocity in the membrane is c = T /σ m .

The boundary condition is y m (a, t) = 0 (membrane clamped at rim). The initial condition representing a stroke with a mallet is ∂y ∂t (r, 0) = -V 0 I ε δ 0 , with amplitude V 0 m s , I ε (r) a smooth function that varies between a value 1 for r = 0 and 0 for r ≥ ε, and δ 0 the Dirac delta function at t = 0.

2) Cavity: Hypotheses (H2a)-(H2b) correspond to a constant value of pV γ , with p the total pressure, V the total volume and γ the heat capacity ratio. We differentiate this relation and approximate the differentials by differences, using acoustical quantities p c and V c that are the deviation values of the instantaneous pressure and volume with respect to the values p 0 and V 0 at rest,

p c = p -p 0 V c = V -V 0 .
One obtains the approximate relation (at order 1)

p c /p 0 = -γV c /V 0 . (2)
The two contributions to the volume change V c are hatched in Figure 1. There is a volume change due to the transverse movement of the upper membrane and a volume change due to the position of the loudspeaker membrane, which is modeled by a plane piston of vertical position l and effective surface S l . Therefore,

V c (t) = S y m (r, t)dS(r) -S l l (t), (3)
with dS(r) = 2πr dr, r ∈ [0, a]. Substituting (3) in ( 2) yields

p c (t) = -γ p 0 V 0 S y m (r, t)dS(r) -S l l (t) . (4)
3) Loudspeaker: The loudspeaker is modeled as a rigid disc of area S l (H3a). The force balance on this disc yields

m l ¨ l (t) + c l ˙ l (t) + k l l (t) = S l p c (t) + Bli l (t). (5)
The left-hand side represents the forces of a standard massspring-damper system. The pressure in the cavity applies a force S l p c (t) on the loudspeaker. The controlled current i l (t) induces a Laplace force Bli l (t), with B the magnetic field [T] and l [m] the length of the loudspeaker coil inside the magnetic field.

D. Spatial eigenmodes

Substituting ( 4) in ( 1) and setting µ = 0 and l (t) = 0 yields a PDE for the movement of an undamped membrane coupled to the cavity without the controlled loudspeaker:

c 2 ∆y m (r, t) - ∂ 2 y m (r, t) ∂t 2 = γp 0 σ m V 0 S y m (r, t)dS(r). (6)
In [START_REF] Morse | Vibration and sound[END_REF], it is shown that for [START_REF] Meurisse | Contrôle actif appliqué aux instruments de musique à vent[END_REF], the displacement y m (r, t) can be written as y m (r, t) = +∞ n=1 ϕ n (r)z n (t) with axisymmetric (under Hypothesis (H1b)) eigenfunctions

ϕ n (r) = w n [J 0 (λ n r) -J 0 (λ n a)] , (7)
where J 0 (•) is a Bessel function of the first kind and w n an arbitrary weight coefficient. It can be shown (see [START_REF] Morse | Vibration and sound[END_REF]) that the wave number λ n satisfies the implicit condition

λ 2 n a 2 J 0 (λ n a) = -πa 4 γp0 T V0 J 2 (λ n a) -BJ 2 (λ n a), (8)
where the dimensionless group B can be interpreted as the ratio of the restoring force applied by the air in the cavity to the membrane tension [START_REF] Kinsler | Fundamentals of acoustics[END_REF].

III. PORT-HAMILTONIAN FORMULATION

In the Port-Hamiltonian System (PHS, [START_REF] Maschke | Port-controlled Hamiltonian systems: modelling origins and systemtheoretic properties[END_REF], [START_REF] Van Der Schaft | Port-Hamiltonian systems: an introductory survey[END_REF]) framework, a system is modeled as consisting of three types of elementary components: energy-storing components, energydissipating components and external components. A PHS model satisfies the power balance. Within this framework, stable simulation [START_REF] Falaize | Passive guaranteed simulation of analog audio circuits: A port-Hamiltonian approach[END_REF] and control [START_REF] Ortega | Interconnection and damping assignment passivity-based control of portcontrolled Hamiltonian systems[END_REF] methods can be developed.

In this section, we recast the system model ( 1)-( 4)-( 5) as a PHS.

A. Definitions

We use following definitions for (in)finite-dimensional PHSs.

Definition 1 (Finite-dimensional PHS): A classical formulation of a PHS of dimension n ≥ 1 and state

x ∈ R n is ẋ = (J -R)∇H(x) + Gu y = G ∇H(x) (9)
where

• u ∈ R q is the input of dimension q ≥ 1,
• y ∈ R q is the output associated with the input (in general different from the measured output y m ), • H(x) ≥ 0 is the Hamiltonian function that represents the total stored energy E of the system as function of its state,

• J = -J is an skew-symmetric (n × n) matrix rep-
resenting the energy exchanges between energy-storing components,

• R = R 0 is a symmetric (n × n) matrix accounting for the dissipation, • G is an (n × q) matrix.
Definition 2 (Infinite-dimensional PHS): An infinite-dimensional PHS is defined as [START_REF] Curtain | An introduction to infinite-dimensional linear systems theory[END_REF], [START_REF] Hélie | Nonlinear damping models for linear conservative mechanical systems with preserved eigenspaces: a port-Hamiltonian formulation[END_REF]

∂ t x = (J -R)δ x H(x) + Gu y = G * δ x H(x) (10)
where

• X is the energy state space,

• u and y are the input and its associated output,

• a scalar product •, • X and norm • X are defined,

• the Hamiltonian function is obtained as

H(x) 1 2 x 2 X ,
• the operator δ x is the variational derivative [START_REF] Villegas | A port-Hamiltonian approach to distributed parameter systems[END_REF],

• the operator J is formally skew-symmetric and the operator R is non-negative symmetric, w.r.t. the scalar product, • G is an operator and G * its adjoint operator.

B. PHS model 1) State space:

Let Ω [0, a]. We define the state as

x = π m (r, t) ε m (r, t) V c (t) π l (t) l (t) ,
whose components are respectively the areal momentum of the membrane (= σ m ẏm ) and its strain (= ∇y m ), the incremental volume of the cavity, the loudspeaker membrane's momentum (= m l ˙ l ) and its vertical position.

The state space is the Hilbert space

X = H 1 0 (Ω, R ≥0 ) × H 1 (Ω, R ≥0 ) × R 3
, endowed with a scalar product x, y X that is defined below according to the total energy of the system.

2) Energy: The total energy of the system is the sum of five contributions.

a) The kinetic energy of the membrane:

1 2 S π 2 m σm dS(r)
b) The elastic energy of the membrane: d 3) Scalar product and Hamiltonian function: We define the scalar product

x, y X S x 1 (r) x 2 (r) L m y 1 (r) y 2 (r) dS(r) + x 3 x 4 x 5 L c&l x 3 x 4 x 5 with L m = diag 1 σ m , T 0, L c&l = diag E 0 V 2 0 , 1 m l , k l 0.
With this scalar product, the Hamiltonian function

H(x) = 1 2
x, x X corresponds to the total energy of the system. 4) Derivative of the Hamiltonian: We are dealing with a mixed PHS, i.e. a combination of a finite-dimensional and an infinite-dimensional PHS. We define a derivative D x H(x) that contains both variational derivatives, for the infinitedimensional state variables, and gradients, for the finitedimensional state variables:

D x H(x) =   δ xm H(x) ∇ x c&l H(x)   = Lx.

Note that

Lx = πm σm T ε E0 V 2 0 V c π l m l k l l ,
respectively the membrane velocity, the radial projection of the membrane's tension force (for small displacements), the cavity energy density, the loudspeaker velocity and the elastic force on the loudspeaker.

5) PHS:

The first equation of the PHS model is (cf. Equation [START_REF] Tiwari | Effects of air loading on the acoustics of an Indian musical drum[END_REF] where δ x is replaced by D x ) is shown in Equation [START_REF] Christian | Effects of air loading on timpani membrane vibrations[END_REF] on the next page. Here, 1 represents the identity operator.

We note that the operator J is formally skew-symmetric 2 , and with domain D(J ) = H 0 (grad) × H(div) × R 3 where grad = ∂ ∂r and div = ∂ ∂r + 1 r in polar coordinates. The output equation of the PHS model corresponds to y = G Lx = Bl ˙ l u l , the electromotive force. Note that the product between input i l and output u l represents the electrical power supplied to the loudspeaker of the system. This output equation of the PHS will not be considered in the remainder of the paper.

∂ ∂t          πm εm Vc π l l          =                   0 div -1 0 0 grad 0 0 0 0 S • dS 0 0 -S l 0 0 0 S l 0 -1 0 0 0 1 0          -          µ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 c l 0 0 0 0 0 0                             πm σm T ε E 0 V 2 0 Vc π l m l k l l           +          0 0 0 Bl 0          i l (11) ∂ ∂ t          πm εm Ṽc πl ˜ l          =                   0 ∂ ∂ r + 1 r -α1 0 0 ∂ ∂ r 0 0 0 0 α S •d S 0 0 -β1 0 0 0 β1 0 -β2 0 0 0 β2 0          -          ξ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ζ 0 0 0 0 0 0                           πm εm Ṽc πl ˜ l          +          0 0 0 1 0          ĩl ( 12
)
d d t          z1 ∂t z1 Ṽc πl ˜ l          =                     0 1 κ 2 0 0 0 -1 κ 2 0 √ πα κBχ 1 λ1sign J0( λ1) 0 0 0 - √ πα κBχ 1 λ1sign J0( λ1) 0 -β1 0 0 0 β1 0 -β2 0 0 0 β2 0           -          0 0 0 0 0 0 ξ κ 2 0 0 0 0 0 0 0 0 0 0 0 ζ 0 0 0 0 0 0                             κ 2 λ2 1 1 -1 B χ -2 1 z1 κ 2 ∂t z1 Ṽc πl ˜ l          +          0 0 0 1 0          ĩl (13)
πm= πm a -1 √ σmγp 0 V 0 εm= εm a -1 √ γp 0 V 0 T -1 Ṽc= Vc V 0 πl = π l √ m l γp 0 V 0 ˜ l = l k -1 l γp 0 V 0 H= H γp 0 V 0 ĩl = i l c(aBl) -1 √ m l γp 0 V 0 ζ= ac l cm l ξ= µa cσm α= a 2 c γp 0 σmV 0 β 1 = aS l c γp 0 m l V 0 β 2 = a c k l m l

C. Dimensionless PHS model

We nondimensionalize the PHS model in order to simplify subsequent calculations, using the transformations shown in Table I.

We thus obtain the dimensionless PHS shown in Eq. ( 12) above with d S = 2πr dr, r ∈ Ω = [0, 1].

IV. MODEL ORDER REDUCTION A. Modal projection

We now project the infinite-dimensional part of the (dimensionless) PHS ( 12) on a spatial basis and truncate at N modes 3 , in order to obtain a finite-dimensional approximative (dimensionless) PHS.

The position of the membrane ỹm is thus approximated as a finite sum ỹm r, t

≈ N n=1 φn (r)z n ( t)
where φn (r) are the spatial eigenfunctions of the membrane coupled to the cavity (7) that form an orthonormal basis 4and zn ( t) are the corresponding temporal evolutions.

(Nondimensionalization was done using ỹm = ym a , φn = ϕ n and zn = 1 a z n .)

B. Port-Hamiltonian approximated model

The result of this projection is shown in ( 13) above for the case of N = 1 mode. (In this formulation, κ a T (γp 0 V 0 ) -1 and χ 1 λ4

1 4B 2 - λ2 1 2B + 1 4 + 2 B .

V. PASSIVE FINITE-TIME CONTROL OF THE UNIMODAL MODEL

A. Control goal and design of the method

Our control goal is to modify the frequency of the first axisymmetric vibration mode (z 1 ) of the tom-tom membrane. To this means, we develop a control loop that controls the loudspeaker current ĩl based on a measurement of the pressure p c inside the cavity and z1 , the temporal evolution of the first mode 5 . By changing the vertical position of the loudspeaker, the volume and thus pressure inside the cavity are changed. This pressure acts on the tom-tom membrane and influences its vibrational behavior.

The synthesis of this controller for the loudspeaker current is based on the projected SHP model for N = 1 mode (13) and consists of 3 steps.

1) The membrane is controlled by a pole placement, the cavity pressure p c being viewed as an intermediate control input.

2) The cavity couples the three subsystems: the temporal evolution z1 of the first membrane mode, the pressure in the cavity p c and the loudspeaker position ˜ l are related.

3) The loudspeaker membrane position is controlled in finite time by the current ĩl in order to achieve the desired pressure p * c ensuring the pole placement of the upper membrane.

B. Control of the membrane by pole placement

The first two equations of the PHS (13) represent the vibration of the first axisymmetric mode of the circular membrane subjected to a distributed pressure pc ( p c /p 0 = -γ Ṽc ). If we are able to control this pressure by a pole placement, i.e.,

p * c t = -k 1 z1 t -k 2 ∂ t z1 t , (14)
with k 1 and k 2 parameters that have been chosen accordingly, we obtain a desired frequency and damping of the first mode. This is possible because the subsystem of the first mode of the membrane controlled by the pressure satisfies the Kalman controllability condition (see e.g. [START_REF] Andréa-Novel | Control theory for engineers -A primer[END_REF]).

C. Coupling between the three subsystems

The third equation of the PHS [START_REF] Krstic | Boundary control of PDEs: A course on backstepping designs[END_REF] provides a relation between the temporal evolution z1 of the first membrane mode, the pressure in the cavity pc and the loudspeaker position ˜ l . It corresponds to

1 γ pc t = √ πακ Bχ1 λ1 sign J 0 λ1 z1 t + β1 β2 ˜ l t . ( 15
)
If we substitute an amplitude measurement zm 1 t and a pressure measurement pm c t in this unimodal relation, we obtain an (implicit) measurement for the loudspeaker position ˜ m l t . If we use a pressure reference p * c t corresponding to the pole placement ( 14) instead, we obtain a reference for the loudspeaker position ˜ * l t .

D. Finite-time control of the loudspeaker

We recall some elements before stating the obtained control law.

1) Definition and finite-time control of a double integrator: A system controlled in finite time will reach an equilibrium point in a finite time. A more formal definition of finite-time stability is given below. Definition 3 (Finite-time stability [START_REF] Bhat | Finite-time stability of continuous autonomous systems[END_REF]): Let ẋ = F (x) represent a closed-loop system, with F (x) continuous and F (0) = 0. Let Φ t (x 0 ) be the time evolution of the state for a given initial state x 0 .

The origin is a finite-time stable equilibrium if there exists an open neighborhood U ⊂ R n of the origin, where following statements hold:

1) Finite-time convergence. There exists a settling-time function T : U \{0} → R ≥0 such that for each x 0 ∈ U \{0}, the evolution Φ t (x 0 ) is defined and unique on t ∈ [0, T (x 0 )[ and lim t→T (x0) Φ t (x 0 ) = 0. 2) Lyapunov stability. There exists a monotonically increasing function δ(•), δ(0) = 0, such that for each x 0 ∈ U , Φ t (x 0 ) ≤ δ( x 0 ) for each t ≥ 0. Furthermore, if U = R n , the origin is globally finite-time stable.

Lemma 1 (Finite-time control of a double integrator [START_REF] Bernuau | Robust finite-time output feedback stabilisation of the double integrator[END_REF]): Consider the double integrator ż1 = z 2 , ż2 = v. The origin is a finite-time stable equilibrium point of this system when it is controlled by the

input v = -k 1 z 1 ξ 2-ξ -k 2 z 2 ξ , with k 1 , k 2 > 0, ξ ∈ ]0, 1[ and x ξ sgn(x)|x| ξ .
2) Finite-time control law for the loudspeaker current: The finite-time control law of Lemma 1 was applied to a mass-spring-damper system in [START_REF] Wijnand | Contrôle des vibrations d'un oscillateur passif : stabilisation en temps fini et par remodelage d'énergie[END_REF]. A passive version of the controller was proposed, in the sense that it can be reinterpreted as a PHS and as such provides robustness against bad parameter estimation. In [START_REF] Lebrun | Electroacoustic absorbers based on finite-time control of loudspeakers: a numerical investigation[END_REF], these considerations were used to formulate a passive law that controls the position ˜ l t of a loudspeaker to a reference trajectory ˜ * l t in finite time. The last equation of the model [START_REF] Krstic | Boundary control of PDEs: A course on backstepping designs[END_REF] represents the (dimensionless) loudspeaker ODE [START_REF] Benacchio | Contrôle actif modal appliqué aux instruments de musique à cordes[END_REF] and can be rewritten as

1 β2 ∂ 2 t ˜ l + ζ β2 ∂ t ˜ l + β 2 ˜ l = -β1 γ pc + ĩl
. By identification, one finds that the corresponding finite-time passive law is

ĩl = β 1 γ pm c +β2 ˜ m l + ζ β 2 ∂t ˜ m l -k3 ˜ m l -˜ * l -k4 ∂t ˜ m l -∂t ˜ * l -1 β 2 k1 ˜ m l -˜ * l ᾱ 2-ᾱ + k2 ∂t ˜ m l -∂t ˜ * l ᾱ , (16)
with control parameters k1 , k2 > 0, ᾱ ∈ ]0, 1[, k3 > β 2 , k4 > ζ β2 . We substitute the signals ˜ m l and ˜ * l (and their time derivatives) as explained below Eq. ( 15) 6 .

Stability: The global system of the membrane and loudspeaker can be viewed in a singular perturbation way. More precisely, the membrane controlled by the pressure can be interpreted as a slow subsystem, exponentially stabilized, whereas the loudspeaker can be interpreted as a fast subsystem, controlled in finite-time through the current. The whole system is then exponentially stabilized according to Tikhonov's theorem (see e.g. [START_REF] De Wit | Theory of Robot Control[END_REF]Thm. A.11]).

VI. SIMULATION

A. Energy-preserving discretization

We simulate the controlled PHS model (13) using a method [START_REF] Hélie | Systèmes Hamiltoniens à Ports avec approche par composants pour la simulation à passivité garantie de problèmes conservatifs et dissipatifs[END_REF], [START_REF] Falaize | Passive guaranteed simulation of analog audio circuits: A port-Hamiltonian approach[END_REF] that ensures stability by preserving the power balance after time discretization by jointly replacing the time derivative of the state by a finite difference and the gradient of the Hamiltonian by a so-called discrete gradient.

B. Simulation

We simulate the unimodal PHS tom-tom model [START_REF] Krstic | Boundary control of PDEs: A course on backstepping designs[END_REF] for which the control law ( 16) is activated at t = 15.

Model parameters: initial state 0 -1 0 0 0 , B = The resulting time evolution of the state is shown in Figure 2 and the periodogram of z1 for the (un)controlled cases is shown in Figure 3. The simulation confirms that the controller is able to modify the frequency of z1 .

Note that the membrane movement exhibits a second, lower frequency in the non-controlled case, due to the coupling with the cavity.

VII. CONCLUSION AND FUTURE WORKS

We have designed a controller for a unimodal model of the membrane. The next step is to extend this result to the multimodal case, considering spill-over, mode interaction and stability of the controller. Furthermore, we wish to implement the controller on a testbench that is under construction.

Fig. 1 .

 1 Fig. 1. Geometry of the system

1 2 S

 2 T ε 2 m dS(r) c) The compression energy of the air in the cavity: Starting from E air = Vc 0 (-p c ) dV c and using the relation (2), one obtains E air = γp0V0

5 .

 5 54, λ1 = 3.08, κ = 0.171, γ = 1.4, α = 1.33, β 1 = 0.768, β 2 = 0.353, ζ = .0923, ξ = 0.1. Control parameters: k 1 = 0.1, k 2 = 0.0005, ᾱ = 0.8, k1 = 10, k2 = 10, k3 = 2β 2 , k4 = 2ζ/β 2 . Simulation parameter: δ t = 0.01.

Fig. 2 .Fig. 3 .

 23 Fig. 2. Simulation of the unimodal tom-tom model: time evolution of the state

We do not consider modes with nodal lines here, for which there is no influence of the cavity.

This can be shown by integration by parts.

Note that the lowest-order modes are most important for the produced sound[START_REF] Kinsler | Fundamentals of acoustics[END_REF]. Furthermore, these modes are influenced the most by the coupling with the cavity[START_REF] Morse | Vibration and sound[END_REF].

After definition of a scalar product φi , φj P S φi φj d S(r) used for the projection, one can prove orthogonality and normalization (given a suited choice of the weights w i ).

In a practical implementation, z1 is deduced from a measurement of the amplitude of the upper membrane, obtained through a laser or piezoelectric sensor.

Note that the practical implementation of this control law based on the measurement of pc and z1 involves additional numerical derivation.