Membrane flutter induced by radiation of surface gravity waves on a uniform flow

Abstract We consider the stability of an elastic membrane on the bottom of a uniform horizontal flow of an inviscid and incompressible fluid of finite depth with free surface. The membrane is simply supported at the leading and the trailing edges which attach it to the two parts of the horizontal rigid floor. The membrane has an infinite span in the direction perpendicular to the direction of the flow and a finite length in the direction of the flow. For the membrane of infinite length we derive a full dispersion relation that is valid for arbitrary depth of the fluid layer and find conditions for the flutter of the membrane due to emission of surface gravity waves. We describe this radiation-induced instability by means of the perturbation theory of the roots of the dispersion relation and the concept of negative energy waves and discuss its relation to the anomalous Doppler effect.


Introduction
Flutter of membranes is a classical subject for at least seven decades.Membranes submerged in a compressible gas flow and their flutter at supersonic speeds have been considered already in the works by Miles (1947), Goland and Luke (1954), Benjamin (1963), and Bolotin (1963).Bolotin (1963), Spriggs, Messiter and Anderson (1969), and Kornecki, Dowell and O'Brien (1976) addressed a problem of the socalled membrane flutter paradox on the relation of stability criteria for an elastic plate to that for a membrane.Absolute and convective hydroelastic instabilities of slender elastic structures submerged in a uniform flow were discussed by Triantafyllou (1992).Recent works on the membrane flutter are motivated by such diverse applications as stability of membrane roofs in civil engineering (Sygulski 2007), flutter of traveling paper webs (Banichuk, Jeronen, Neittaanmäki and Tuovinen 2010;Banichuk, Barsuk, Jeronen, Tuovinen and Neittaanmäki 2019), aerodynamics of sails and membrane wings of natural flyers (Newman and Paidoussis 1991;Tiomkin and Raveh 2017), as well as design of piezoaeroelastic systems for energy harvesting (Mavroyiakoumou and Alben 2020).
Surface gravity waves on a motionless fluid of finite depth is a classical subject as well, going back to the seminal studies of Russell and Kelvin (Carusotto and Rousseaux 2013).Numerous generalizations are known taking into account, for instance, a uniform or a shear flow and surface tension (Maissa, Rousseaux and Stepanyants 2016), flexible bottom (Mohapatra and Sahoo 2011) or a flexible plate resting on a free surface (Das, Sahoo and Meylan 2018a,b,c;Bochkarev, Lekomtsev and Matveenko 2016;Schulkes, Hosking and Sneyd 1987).The latter setting has a straightforward motivation in dynamics of sea ice and a less obvious application in analogue gravity experiments (Carusotto and Rousseaux 2013;Barcelo, Liberati and Visser 2011;Weinfurtner, Tedford, Penrice, Unruh and Lawrence 2011).Recent work (Robertson and Rousseaux 2018) discusses effects of viscous dissipation of surface gravity waves to the analogue gravity.
Remarkably, another phenomenon that is being analysed from the analogue gravity perspective is super-radiance (Carusotto and Rousseaux 2013;Barcelo, Liberati and Visser 2011;Brito, Cardoso and Pani 2015) and its particular form, discovered by Ginzburg and Frank (Ginzburg and Frank 1947;Ginzburg 1996), known as the anomalous Doppler effect (ADE) (Bekenstein and Schiffer 1998;Nezlin 1976;Nemtsov and Eidman 1987).In electrodynamics, the ADE manifests itself when an electrically neutral overall particle, endowed with an internal structure, becomes excited and emits a photon during its uniform but superluminal motion through a medium, even if it started the motion in its ground state; the energy source is the bulk motion of the particle (Bekenstein and Schiffer 1998).
Anomalous Doppler effect in hydrodynamics was demonstrated for a mechanical oscillator with one degree of freedom, moving parallel to the border between two incompressible fluids of different densities (Gaponov-Grekhov, Dolina and Ostrovskii 1983).It was shown that the oscillator becomes excited due to radiation of internal gravity waves if it moves sufficiently fast.In (Abramovich, Mareev and Nemtsov 1986) the ADE for such an oscillator was demonstrated due to radiation of surface gravity waves in a layer of an incompressible fluid.Nemtsov (1986) was the first who considered flutter of an elastic membrane resting on the bottom of a uniform horizontal flow of an inviscid and incompressible fluid as an anomalous Doppler effect due to emission of long surface gravity waves.In the shallow water approximation, he investigated both the case of a membrane that spreads infinitely far in both horizontal directions and the case when the width of the membrane in the direction of the flow (or the chord length) is finite whereas the span in the perpendicular direction is infinite.Nevertheless, the case of the flow of arbitrary depth has not been studied in (Nemtsov 1986) as well as no numerical computation supporting the asymptotical results has been performed.Another issue that has not been addressed in (Nemtsov 1986) is the relation of stability domains for the membrane of the finite width to that for the membrane of the infinite width.Vedeneev (2004) studied flutter of an elastic plate of finite and infinite widths on the bottom of a uniform horizontal flow of a compressible gas occupying the upper semispace.He performed analysis of the relation of stability conditions for the finite plate with that for the infinite plate using the method of global stability analysis by Kulikovskii (Doaré and de Langre 2006;Vedeneev 2016).A single-mode high frequency flutter due to a negative aerodynamic damping and a binary flutter have been identified in (Vedeneev 2016).However, no connection has been made to the ADE and the concept of negative energy waves.
In the present work we reconsider the setting of Nemtsov in order to address the finite height of the fluid layer, find flutter domains in the parameter space, analyze them using perturbation of multiple roots of the dispersion relation and investigate the flutter onset for the membrane of infinite chord length.We will explain the radiative instabilities via the interaction of positive and negative energy waves using an explicit expression for the averaged total energy derived rigorously from physical considerations and relate them to the anomalous Doppler effect.

Physical system
In a Cartesian coordinate system OXY Z consider an inextensible elastic rectangular membrane strip of constant thickness h and density ρ m , of infinite span in the Y -direction, held at Z = 0 at the leading edge (X = 0) and at the trailing edge (X = L) by simple supports, Fig. 1.
The membrane is initially still and flat, immersed in a layer of an inviscid, incompressible fluid of constant density ρ with free surface at the height Z = H.The two-dimensional flow in the layer is supposed to be irrotational and moving steadily with velocity v in the positive X-direction.The bottom of the fluid layer at Z = 0 is supposed to be rigid and flat for X ∈ (−∞, 0] ∪ [L, +∞).Nemtsov (1986) assumed that vacuum is below the membrane.In the present study we prefer to consider that a motionless incompressible medium of the same density ρ is present below the membrane with the pressure that is the same as the unperturbed pressure of the fluid (Vedeneev 2004(Vedeneev , 2016)).
Assuming small vertical displacement of the membrane w(X, t), where t is time, a constant tension T along the membrane profile, and neglecting viscous forces, we write the dimensional membrane dynamic equation as (Tiomkin and Raveh 2017) where ∆P (X, t) is the pressure difference across the interface Z = 0.The simplysupported boundary conditions for the membrane yield In general, to recover the pressure P (X, Z, t) of the fluid we write the Euler equation for the vorticity-free flow (Carusotto and Rousseaux 2013;Maissa, Rousseaux and Stepanyants 2016) with v = ve X + ∇ϕ, where ϕ(X, Z, t) is the potential of the fluid, e X is the unit vector in the X-direction, and g stands for the gravity acceleration.This yields the integral of Bernoulli P ρ (2.4) The incompressibility condition takes the form From the equation (2.4) it follows that in the case when a motionless medium of density ρ is present below the membrane with its pressure equal to the unperturbed pressure of the fluid above the membrane, the linear in ϕ expression for the pressure difference, ∆P (X, t), is (2.6) For the sake of completeness, we present also the analogous expression for the pressure difference for the case when vacuum is below the membrane (Nemtsov 1986) (2.7) Impermeability of the rigid bottom implies the condition (2.8) The prescription of normal velocity at the boundaries of moving surfaces allows us to express the kinematic condition for the membrane and to specify the same condition at the free surface where u(X, t) is the free surface elevation and n is the vector of the outward normal to a surface.This implies that the projection of the vector ∇ϕ to the normal will coincide with the pozitive z-direction for the free surface and have the opposite direction for the membrane, see Fig. 1.Using the Bernoulli integral (2.4) at the free surface where P = 0 and retaining only linear in ϕ terms, we find gu = −(∂ t + v∂ X )ϕ. (2.11) Expressing u from (2.11) and substituting it into (2.10)we obtain the boundary condition at the free surface of the liquid that reads (2.12) 2.2.Dimensionless mathematical model Let us choose the height of the fluid layer, H, as a length scale, and ω −1 0 , where ω 0 = g/H, as a time scale.Then, we can introduce the dimensionless time and coordinates (2.14) the dimensionless parameters of the added mass ratio (Minami 1998) and membrane chord length and the two Mach numbers (Vedeneev 2004(Vedeneev , 2016) where c 2 = T /(ρ m h) is the squared speed of propagation of elastic waves in the membrane and √ gH is the speed of propagation of long surface gravity waves in the shallow water approximation.The added mass ratio α is the ratio of the fluid and solid mass contained, respectively, in the volume delimited by the dashed lines in Fig. 1 and in the membrane (Minami 1998).In Fig. 1, Ω denotes the fluid domain and ∂Ω 0 , ∂Ω 1 , ∂Ω 2 stand, respectively, for the free surface, membrane, and solid wall borders.
The dimensionless wave equation (2.1) is therefore (2.17) Supplementing it with the expression (2.6), which in the dimensionless time and coordinates has the form we find The dimensionless boundary condition (2.9) is whereas the boundary condition (2.12) at the free surface in the dimensionless form looks like (2.21) Collecting together the equations (2.19), (2.20), (2.21) and the obvious dimensionless versions of the equations (2.5) and (2.8) and assuming a time dependence in the form of φ, ξ ∼ e −iωτ results in the following dimensionless set of equations and their boundary conditions for the case when a motionless medium is present below the membrane: ] is the impermeability condition for the membrane.For simplicity, we retain the same notation for the membrane displacement and the fluid potential after the separation of time.Therefore, due to the irrotational, incompressible and inviscid character of the fluid, our mathematical model (2.22) consists of the Laplace equation for the fluid potential (2.22a) supplemented with the kinematic conditions for the free surface (2.22b) and the membrane (2.22c).The pressure at the surface of the fluid is also prescribed as a dynamic condition and therefore closes the system of equations for the fluid in this model: the motion of the membrane is described by a nonhomogeneous wave equation (2.22e) with the pressure of the fluid (recovered through the Bernoulli principle) as a source term.The membrane is supposed to be simply supported at its extremities (2.22f).

Dispersion relation for the fluid layer of arbitrary depth
In this paper, we analyse the case when the chord of the membrane is infinite, i.e. when the membrane extends from −∞ to +∞ in the x-direction.Since the motion of the fluid is two-dimensional in the (x, z)-plane and the horizontal extension of the fluid layer is infinite in the x-direction too, we can represent the potential of the fluid φ in the physical space by means of the inverse Fourier transform of the potential φ in the wavenumber space as where κ is the wavenumber and under standard assumption that both φ(x) and φ(κ) are absolutely integrable functions, implying their vanishing at infinity.Assuming that ∂ x φ is also absolutely integrable, which allows us to use twice the property ∂ x φ = iκ φ, we find the Fourier transform of the Laplace equation (2.22a) The general solution to equation (3.4) where A(κ, ω) and B(κ, ω) are yet to be determined from the Fourier-transformed boundary conditions.The boundary condition (2.22c), expressing the impermeability of the membrane at z = 0, takes the form because the outward direction of the normal vector n to the surface of the membrane is opposite to the positive z-direction, see Fig. 1.The Fourier transform of (3.5) reads Similarly transforming the boundary condition (2.22b) at the free surface we find (3.9) Solving equations (3.8) and (3.10) simultaneously with respect to A and B, we obtain The Fourier transform of the non-homogeneous wave equation (2.22e) for the membrane displacement evaluated at z = 0 reads Inserting the expression (3.4) for φ with the coefficients (3.11) into (3.12),discarding ξ in the result and introducing the new parameters, the phase velocity σ = ω/κ and the coupling parameter β = α/κ 2 , we obtain the following dispersion equation in the case where a medium with constant pressure is present below the membrane It is instructive to show another way of deriving the dispersion equation (3.13).For this, we notice that (3.6) and (3.7) allow us to express ξ by means of ∂ z φ.Using the result in (3.12), we can obtain a boundary condition for φ(z) at z = 0.This new boundary condition together with the boundary condition (3.9) and the equation (3.3) produce a closed-form boundary value problem for the Laplace equation with the Robin boundary conditions: Substituting the general solution (3.4) into the boundary conditions of the problem (3.14) results in the system of two linear equations with respect to A and B This system can be written in the matrix form as where the 2 × 2 matrices involved are Computing the determinant of the matrix polynomial we arrive at the dispersion equation which, with the notation σ = ω/κ and β = α/κ 2 , transforms exactly to (3.13).
For the sake of completeness we present also the dispersion relation for the system with vacuum below the membrane In the shallow water approximation corresponding to the limit κ → 0, the expression (3.19) reduces to .20) which is nothing else but the shallow water dispersion relation derived in (Nemtsov 1986).
In order to get the dispersion relation (3.19), one must take the pressure difference (2.7), make it nondimensional and use in the expression (2.17) which then reads as After separation of time it reduces to the analogue of the boundary condition (2.22e) Inserting the expression (3.4) for φ with the coefficients (3.11) into (3.21)results after familiar algebraic manipulations in the dispersion relation (3.19).

Analysis of the dispersion equation
In the absence of coupling between the free surface and the membrane, i.e. for β = 0, both the dispersion relation (3.13) and the dispersion relation (3.19) reduce to (ω − κF r) 2 = κ tanh κ after taking into account that σ = ω/κ and that M as it is defined in (2.16) can also be interpreted as the Froude number, F r.
The roots of the decoupled dispersion equation (3.22) are real If we consider the roots (3.23) as functions of the fluid Mach number, M , we find that σ ± 1 are two horizontal straight lines and σ ± 2 are two straight lines with the slope equal to 1, see the leftmost panels in Fig. 2. One can see that at β = 0 the root branches intersect at four points forming the double roots σ 0 = M w at and the double roots 24) and (3.23) is the condition of 'phase synchronism' for the case of arbitrary height of the fluid layer that extends the corresponding result obtained in (Nemtsov 1986) in the shallow water limit, κ → 0.
With the increase in β the roots ±σ 0 situated at M = ±M − 0 split into simple real ones and this splitting is accompanied by unfolding the crossings into avoided crossings, Fig. 2. Quite in contrast, the roots ±σ 0 situated at M = ±M + 0 split into complex-conjugate pairs that form bubbles of instability at moderate values of β that open up with the increase in β to originate disclosed complex branches, as is seen in the rightmost panels of Fig. 2.
Let us re-write the dispersion relation (3.13) as follows (3.25)Then, we can apply to it the perturbation theory derived in the Appendix A.
Consider the double root σ 0 at M = M + 0 and β = β 0 = 0. Adapt the approximate equation (A 16) to our model where ∆σ = σ − σ 0 , ∆M = M − M + 0 , and ∆β = β.Calculating the partial derivatives at σ = σ 0 , M = M + 0 , and β = β 0 = 0, we find (3.27) (3.28) For any β > 0 the crossing of the real roots σ at M = M + 0 unfolds into two hyperbolic branches of the real roots that are connected to the "bubble" of complex eigenvalues with the real parts Re σ = and with the imaginary parts that form an ellipse in (M, Im σ)plane  crossing point (3.31) The bubble of instability (3.30) corresponds to the inner points of the instability domain bounded by (3.31).
Using the same methodology to approximate the avoided crossing close to M = M − 0 , σ = σ 0 and β = β 0 by equation (3.26), we obtain Separating real and imaginary parts of σ in (3.32) similarly to how it has been done in the previous case, one can see that the bubble of instability does not originate for β > 0 in the unfolding of the crossing at M = M − 0 , see Fig. 3.In Fig. 4 we show that the exact neutral stability boundaries obtained from equating the discriminant of the fourth-order polynomial (3.25) in σ to zero and their approximation (3.31) calculated at the crossing point at M = M + 0 are in a very good agreement.It is instructive to change the point of view and to look at the critical values of parameters as functions of the Mach number M w of elastic waves in the membrane.In Fig. 5 we present stability maps of the dispersion equation (3.25) given by its discriminant in the (M w , β)-plane for the fixed value of M = M 0 = 2 and increasing values of κ.We see that for all κ the instability is possible only in the interval |M w | < M 0 = 2, which agrees with Fig. 4. For β = 0, the instability domain touches the M w -axis at the points M w = M 0 − tanh κ κ and M w = −M 0 + tanh κ κ .In the limit κ → 0, the touching occurs at M w = M 0 − 1 = 1 and M w = −M 0 + 1 = −1.
A qualitative change happens when κ κ 0 where κ 0 > 0 is uniquely determined by M 0 > 0 from the equation For instance, M 0 = 2 yields κ 0 ≈ 0.5218134478.At κ = κ 0 a new, isolated, domain of instability originates that touches the M w -axis at β = 0 and is growing when κ is further increasing, Fig. 5.At some value of κ the two domains touch each other and then form a unified domain.At κ → ∞ the central part of the unified domain dominates over its side parts corresponding to the instability found by Nemtsov in the shallow-water approximation when κ → 0 and the coupling β is weak, Fig. 5.
To understand the origin of the new instability, we plot the real and imaginary values of σ as functions of M w in Fig. 6 for a given M = M 0 = 2.The central panel of Fig. 6 corresponding to β = 0 and κ = κ 0 shows four straight lines intersecting at five points, including the origin.The upper horizontal line corresponds to the fast surface gravity wave with σ = M 0 + tanh κ0 κ0 ≈ 3 whereas the lower horizontal line to the slow surface gravity wave (Nemtsov 1986) with σ = M 0 − tanh κ0 κ0 ≈ 1.The two inclined lines correspond to the forward and backward elastic waves in the membrane with σ = ±M w .When β and κ depart respectively from zero and κ 0 , all the five crossings unfold either into avoided crossings (as the elastic and fast surface gravity waves) or into bubbles of instability (as the elastic and slow surface gravity waves) resulting in the high-frequency flutter due to radiation of long surface gravity waves.For β > 0 the crossing at the origin transforms into an avoided crossing at κ < κ 0 or into the bubble of instability at κ > κ 0 , which yields low-frequency flutter at short wavelengths κ.
The left panel in Fig. 7 allows us to track the evolution of the flutter domains as κ varies from zero to infinity at M 0 = 2 and β = 0.03.Nemtsov's radiation-induced flutter domain is the widest in the shallow water limit and evolves along the curves (shown as black solid lines in Fig. 7(left)) to which the Nemtsov domains degenerate at β = 0. Note that the Nemtsov flutter domain is perfectly approximated by the formula (3.31) obtained from the unfolding of the eigenvalue crossing corresponding to the slow surface gravity wave and the elastic wave (dotted lines in Fig. 7(left)).
To understand the central instability domain shown in Fig. 7(left) for a given β we plot it in the (M w , κ, β)-space in the right panel of Fig. 7, given M = M 0 .One can see that the domain is symmetric with respect to the plane M w = 0 and has a pronounced conical singularity at κ = κ 0 determined by the equation (3.33) when β = 0 and M w = 0. Equation (3.33) follows from the discriminant of the dispersion equation (3.25) at β = 0 and M w = 0.The conical singularity of the stability boundary therefore exactly corresponds to the crossing of the eigenvalue curves at the origin in the central panel of Fig. 6.Usually, the conical singularity of the stability boundary is associated with a double semi-simple eigenvalue with two linearly-independent eigenvectors (Kirillov and Seyranian 2004;Kirillov 2009Kirillov , 2010Kirillov , 2013;;Guenther and Kirillov 2006;Kirillov, Guenther and Stefani 2009).
By this reason, we apply the perturbation theory of double eigenvalues presented in Appendix A to the double zero eigenvalue σ = σ 0 = 0 at the crossing shown in the central panel of Fig. 6 and corresponding to the values of parameters β = β 0 = 0, κ = κ 0 , M = M 0 , M w = M w,0 = 0.A natural extension of the approximation formula (A 16) to the case of four parameters β, κ, M w , and M yields Computing the corresponding partial derivatives of the left part of the dispersion equation (3.25), and evaluating them at β = β 0 = 0, κ = κ 0 , M = M 0 , M w = M w,0 = 0, where M 0 and κ 0 are related by the equation (3.33), we find that the only non-zero derivatives are (3.36)Taking this into account in (3.35), we find a simple approximation describing the unfolding of the double zero eigenvalue (3.37) Let us further assume that M = M 0 is fixed.Then the last term in (3.37) vanishes, and the discriminant of the resulting quadratic polynomial in σ produces the equation of a cone with the apex at κ = κ 0 , M w = 0, and β = 0 38) The cone (3.38) is shown in red in the right panel of Fig. 7.With β = 0.03, M 0 = 2 and κ 0 computed by means of the equation (3.33), the approximation (3.38) fits the boundary of the exact instability domain with a remarkable precision, as is evident in Fig. 7(left).
It is easy to see that in the plane M w = 0 the cone (3.38) defines the two lines that approximate the instability domain near κ = κ 0 , see the central panel in Fig. 8.
As soon as M w deviates from zero, the cone (3.38) again provides a very good fit to the actual stability boundary, Fig. 8(right).In the plane the cross-section of the cone (3.38) is described by the two lines that constitute a linear approximation to the stability boundary shown in Fig. 8(left).

Wave energy of the Nemtsov system with the membrane of infinite chord length
Let us use physical considerations to derive the expression for the averaged over the wave period energy of the Nemtsov system with the membrane of infinite chord length, by combining the approaches of the works by Maissa, Rousseaux and Stepanyants (2016) and Schulkes, Hosking and Sneyd (1987).
In the linear wave theory, the energy is a function of the squared wave amplitude (Maissa, Rousseaux and Stepanyants 2016).Therefore, the total energy per surface area of the membrane resulting both from the wave velocity of the structure and the elastic energy due to its tension is (3.42)where Re stands for the real part of the vibration amplitude that is complex-valued because of the assumed plane wave solution Recall that φ(z) is determined by the expression (3.4) with the coefficients (3.11) and η, ξ are respectively displacement amplitudes of the free surface and the membrane.
The energy of the fluid depends on whether we assume a vacuum below the membrane (Nemtsov 1986) or a motionless medium of the same density as the fluid above the membrane and with the pressure equal to the pressure of the unperturbed fluid (Vedeneev 2004(Vedeneev , 2016)).The gravitational potential energy of the free surface is the only term contributing to the total potential energy of the fluid in the latter context.Therefore, The kinetic energy of the flow per unit area is determined by the velocity field u = ∇φ + M e x , where u = v √ gH , that needs to be directly integrated within the limits given by the surface of the membrane and the free surface of the fluid Similarly, with the help of the integration by parts, the Laplace equation (3.3), and the expressions (3.46), we obtain Finally, following (Maissa, Rousseaux and Stepanyants 2016), we evaluate the last integral term in (3.45) with the help of the Lagrange mean value theorem, which is justified by the assumption that η and ξ are infinitesimally small perturbations of the surface boundaries ∂Ω 0 and ∂Ω 1 .Performing this procedure, and then taking into account expressions (3.46), we derive (3.49) Notice that the right hand sides in the expressions (3.47), (3.48), and (3.49) are T -periodic functions of time, where T = 2π/ω.Averaging these expressions over the wave period T according to the rule we deduce the mean kinetic energy of the fluid The term φ∂ z φ in (3.51) is evaluated with the help of the Bernoulli principle (2.11) and the free surface kinematic condition (2.10) at ∂Ω 0 , and the wave equation (2.22e) with the impermeability condition (2.22c) at ∂Ω 1 .This yields, respectively, Substituting expressions (3.52) into (3.51)we obtain the final expression for the mean kinetic energy of the fluid The other energies of the system, after taking into account (3.46) and time-averaging (3.50), become Notice that in the absence of the background flow (M = 0) the system respects the equipartition of energy in accordance with the virial theorem (Landau and Lifschitz 1987), because the flow is irrotational and thus derived from a fluid potential (Schulkes, Hosking and Sneyd 1987).
After summing up all the different terms given by the equations (3.53) and (3.54) we obtain the total averaged energy thus providing an extension to the case when the velocity field contains a background flow (M = 0).A more suitable expression for the mean total energy can be obtained by expressing the different amplitudes of the system in (3.55) in terms of a unique one, for instance, ξ.From the kinematic condition (2.10) on the free surface with the plane wave solution Next, expressing the term κ 2 M 2 w by means of the dispersion relation (3.18) and substituting the result into (3.57)yields a more compact formula for the total energy: (3.58) Notice that the term in the braces in Eq. (3.58) is nothing else but the partial derivative ∂D ∂ω of the dispersion relation (3.18) written in the following equivalent form (3.59) where D m = ω 2 −κ 2 M 2 w and D f = [(ω −κM ) 2 −κ tanh κ] stand for the dispersion relation of, respectively, the free membrane and the free surface flow with a rigid boundary at the bottom.This proves that our total energy per unit area, averaged over the wave period, (3.60) The representation (3.60) can be found, e.g. in (Cairns 1979), and can be derived in the frame of the general Lagrangian variational approach (Whitham 1999;Ostrovskii, Rybak and Tsimring 1986), see also the recent work (Fukumoto, Hirota and Mie 2014) for historical notes and application to stability of vortices.Notice that according to (3.60) the energy vanishes at the points where ω = 0 or ∂D ∂ω = 0, the latter condition corresponding to the existence of multiple roots of the dispersion relation.Correspondingly, the ratio E /ω, which is the averaged wave action A (Zhang, Qin, Davidson, Liu and Xiao 2016), vanishes only at the locations of the multiple eigenvalues, cf.Fig. 2

Discussion
Comparing the eigenvalue plots of Fig. 2 and Fig. 6 with the corresponding to each branch averaged wave energy and wave action that are shown in Fig. 9 and Fig. 10, respectively, we notice that flutter instability is necessarily accompanied with the interaction of waves of opposite sign of energy/action.In contrast to the action, the energy changes sign also at the points where the phase velocity σ changes sign, quite in accordance with (3.60).
Looking now at the roots (3.23) of the decoupled dispersion equation (3.22), we conclude that the elastic waves σ ± 1 = ±M w propagating in the membrane always have positive energy whereas among the surface gravity waves σ ± 2 = M ± tanh κ κ it is the energy of the slow wave σ − 2 that becomes negative for M > 0 as soon as M > tanh κ κ .Therefore at the crossing (3.24) corresponding to M + 0 = M w + tanh κ κ the positiveenergy/action elastic wave meets the slow surface gravity wave that carries negativeenergy/action (Nemtsov 1986).
With β increasing from zero the crossing unfolds causing the eigenvalue branches to merge on the interval bounded by the points where ∂ ω D = 0.At these exceptional points (Kirillov 2013) both the energy and the action change sign, see Fig. 9 and Fig. 10.On the interval the roots are complex and form the bubble of instability, see Fig. 2 and Fig. 6.
Since the fast surface gravity wave carries positive energy, one needs to add energy to the flow in order to excite this wave on the flow.In contrast, in order for the slow surface gravity wave carrying negative energy to build up on the flow, the energy must be extracted from the flow (Nezlin 1976) via some mechanism for dissipation of its energy.
In the Nemtsov problem, such a mechanism is the energy transfer from the slow surface gravity wave to an elastic wave associated with the membrane, which is a stationary medium and therefore has positive energy (Nezlin 1976).One can say that this transferred energy yields flutter of the membrane due to emission of the slow surface gravity wave carrying negative energy.In Fig. 2 and Fig. 6 as well as in Fig. 9 and Fig. 10 we observe that the flutter instability of the membrane occurs only if the velocity of the flow is higher than the phase velocity of the oscillations on the surface of the flow, σ < M , i.e. the flow moves faster than the waves it can excite (Nemtsov and Eidman 1987;Nezlin 1976).The condition ω = M κ or σ = M is known as the Cerenkov condition for emission of radiation by a moving source (Nezlin 1976;Ginzburg and Frank 1947;Ginzburg 1996;Bekenstein and Schiffer 1998;Carusotto and Rousseaux 2013).Being substituted into a dispersion relation, the Cerenkov condition transforms the latter into the expression defining a surface in the space of wave numbers that determines the wake pattern behind the source (Carusotto and Rousseaux 2013;Schulkes, Hosking and Sneyd 1987).For the supercritical velocities M > σ the surface in the space of wave numbers develops a conical singularity known as the Cerenkov cone (Carusotto and Rousseaux 2013;Nemtsov 1986) with the angular aperture 2 arccos σ M .
Anomalous Doppler effect (ADE) is the change in the sign of the field frequency radiated into the Cerenkov cone as compared with the field radiated outside this cone (Nezlin 1976;Nemtsov 1986;Gaponov-Grekhov, Dolina and Ostrovskii 1983;Abramovich, Mareev and Nemtsov 1986;Carusotto and Rousseaux 2013).It is exactly the slow surface gravity wave that satisfies this condition Hence, for the one-dimensional or, more precisely, plane two-dimensional waves, both the negative energy waves and the ADE correspond simply to waves with phase velocity lower than the flow velocity and wave-vector pointing in the same direction as the flow (Ostrovskii, Rybak and Tsimring 1986;Nemtsov 1986).In our case, the radiated slow gravity wave increases the energy of oscillations of the membrane at the expense of the energy of the flow that supports this wave.
Finally, we plot the dispersion curves ω(κ) in Fig. 11, which show that the slow surface gravity wave branch and the membrane branch interact under the line Re(ω) = κM , if κ > 0. Substituting the Cerenkov condition in the dispersion relation (3.18) we reduce it to (M 2 w − M 2 ) tanh κ = 0, thus providing a rationale for the absence of instabilities for M 2 w > M 2 that is evident in all our stability diagrams.

Conclusion
Through the revival of a classical hydrodynamical model performed in this work, we have been able to extend the stability analysis and to enhance knowledge of the underlying physics by making connections with the fundamental concepts such as negative energy waves and the anomalous Doppler effect, supported by advanced mathematical tools.
Our continuation of Nemtsov's investigation of the radiation-induced flutter of a membrane in a uniform flow with the new derivation of the dispersion relation for the fluid layer of arbitrary depth and membrane of infinite chord length has led to a significant improvement in the computation of stability diagrams without any limitation on the range of the system parameters.
The exhaustive stability analysis has been performed using the original perturbation theory of multiple roots of the dispersion relation to compare with the exact stability domains and both computations are proven to be in an excellent agreement.More precisely, the crossings and avoided crossings of the dispersion curves are very well approximated by the simplified expressions for the phase speed of the membrane and fluid modes derived with the perturbation approach.
After computing the discriminant of the full dispersion relation, we have identified a new instability domain arising from a conical singularity in the parameter space that was not present in the early study of Nemtsov.This new domain is associated with a low-frequency flutter for short wavelengths and corresponds to the case when the velocity of propagation of elastic waves in the membrane is much smaller than the velocity of the flow.
Moreover, following the procedures used in previous studies on simplified hydrodynamical systems to calculate the averaged wave energy and after developing the method further to take into account the coupling between the free surface of the flow and the elastic membrane on the bottom, we have obtained an elegant and applicable expression for the total averaged energy.We have verified that in the absence of the background flow, the system respects the equipartition of energy in accordance with the virial theorem, thus confirming that the existence of the negative energy waves can only be possible when the flow is in motion.
We have shown that the formula for the total averaged energy recovered in our work by means of the direct integration of different physical fields is expressed via the derivative of the dispersion relation with respect to the frequency of oscillations and reduces exactly to the form described by Cairns (1979).
The anomalous Doppler effect (ADE) is a direct consequence of the relative motion of an oscillator in a medium and more precisely, it occurs when the internal energy of the system increases due to the emission of negative energy waves (NEW).In our context, while the system is composed of a fluid layer and a membrane, such phenomenon has been proved by Nemtsov to exist only when the conditions of phase synchronism and NEW emission are satisfied.The criterion for the phase synchronism in the system is easily identified in the computations of the dispersion curves as the crossings of the different branches that lead to the onset of positive growth rate and therefore to temporal instability.The latter phenomenon is a natural consequence of the highly excited state of energy that the system transit to due to the dominance of NEW over the waves carrying positive energy.Indeed, NEW emission is known as a process that increases the total energy of an oscillatory system while radiating energy away from the oscillator and, only when this gain in internal energy exceeds the loses from the contribution of positive energy modes, the total energy of the system starts growing in amplitude.Hence, it requires to have waves carrying energy of opposite signs that interact for the instability to develop.
Our expressions for the action and energy of the Nemtsov system demonstrated as expected the collision of waves carrying positive and negative energy as the onset for the radiative instability and the flutter of the membrane.Such phenomenon is well-known in the community of dynamical systems but in this context, it is associated with the emission of NEW in the region of ADE.Hence, in addition to improving the stability

Figure 1 .
Figure 1.An elastic membrane with the chord of length L attached to two rigid walls along its leading (X = 0) and trailing (X = L) edges on the bottom of a fluid layer of depth H moving with the velocity v. Ω is the fluid domain and ∂Ω0, ∂Ω1, ∂Ω2 are respectively the free surface, membrane, and rigid wall borders.

Figure 3 .
Figure 3.For Mw = 1, κ = 1, and β = 0.1 real and imaginary parts of the roots of the dispersion relation (red) (3.13) and (blue, dashed) their approximations by the equations (3.28) and (3.32) near the crossing points that exist at β = 0, M = M ± 0 , σ = σ0.Notice an avoided crossing above the line Re(σ) = M and the bubble of instability below this line.

Figure 4 .
Figure 4. Stability maps of the dispersion equation (3.25) given by its discriminant for (upper left) Mw = 1 and κ = 1, (upper right) Mw = 1 and β = 0.1, (lower left) Mw = 1 and M = 1.6, and (lower right) β = 0.5 and κ = 1.The regions of real phase speed σ are shown in white (stability) and that of the complex σ (temporal instability) in light blue.The red dotted curves is the approximation (3.31).Notice the absence of instabilities for Mw > M in the lower right panel.

Figure 5 .
Figure 5. Stability maps of the dispersion equation (3.25) given by its discriminant for M = M0 = 2 and (from upper left to lower right) κ = 0.5, κ = 0.55, κ = 0.58, κ = 0.8, κ = 1.5, κ = 3.The regions of real phase speed σ are shown in white (stability) and that of the complex σ (temporal instability) in light blue.The black dotted curves is the approximation (3.31) and the solid red ellipse is the conical approximation (3.38).When κ → ∞, the central part of the instability domain approximated by (3.38) dominates over the side parts of the domain.Notice the absence of instabilities for Mw > M0.
.30) see Fig.3.Equating to zero the discriminant of the quadratic in σ equation (3.28), we arrive at the following quadratic approximation to the neutral stability curve at the

Figure 7 .
Figure 7. (Left) Stability map of the dispersion equation (3.25) given by its discriminant for M = M0 = 2 and β = 0.03.The regions of real phase speed σ are shown in white (stability) and that of the complex σ (temporal instability) in light blue.The black dotted curves is the approximation (3.31) and the solid red line is the conical approximation (3.38).When β = 0, the blue instability domains degenerate (central) to the ray κ κ0 ≈ 0.5218134478 and (sides) to the curves (3.34) shown as solid black lines.(Right) Stability boundary with the conical singularity at κ = κ0, β = 0 and Mw = 0, according to (blue, internal surface) the discriminant of the dispersion equation (3.25) and (red, external surface) to the approximation of the cone (3.38).

Figure 8 .
Figure 8.For M = M0 = 2 (left) cross-section of the instability domain with the conical singularity shown in Fig. 7(right) by the plane (3.40).The regions of real phase speed σ are shown in white (stability) and that of the complex σ (temporal instability) in light blue.The red lines crossing at the apex of the cone at κ = κ0 ≈ 0.5218134478 is a linear approximation given by Eq. (3.41).(Centre) Cross-section by the plane Mw = 0 of the instability domain and (red line) its linear approximation (3.39) at the conical point κ = κ0.The black dotted line is given by the Eq.(3.40).(Right) Similar cross-section by the plane Mw = 0.1 where the red curve is the approximation (3.38).
and Fig 9.In the latter figure as well as in Fig. 10 we show several computations of the averaged wave energy and wave action over the fluid Mach number M , and, respectively, the membrane Mach number M w , where ω is calculated with the use of the dispersion relation (3.18).