J Sierra

D Fabre

V Citro

F Giannetti

Acoustic instability prediction of the flow through a circular aperture in a thick plate via an impedance criterion

Keywords: Instabilities, acoustics, jet, impedance

We investigate the mechanisms leading to acoustic whistling for a jet passing through a circular hole in a thick plate connecting two domains. Two generic situations are considered. In the first one, the upstream domain is a closed cavity while the downstream domain is open, leading to a class of conditionally unstable modes. In this case, the instability source lies in the recirculation region within the thickness of the plate, but coupling with a conveniently tuned resonator is needed to select the conditional instability range. In the second situation, the two regions, upstream and downstream of the hole, are considered as open, leading to a class of hydrodynamic modes where instability of the recirculation region is sufficient to generate self-oscillations without the need of any resonator. A matched asymptotic model, valid in the low-Mach limit, is used to derive a global impedance of the system, combining the impedance of the hole and the modelled impedances of the upstream and downstream domains. It is shown that the knowledge of this global impedance along the real ω-axis provides an instability criterion and a prediction of the eigenvalues of the full system. Validations against the solutions of the eigenvalue problem obtained from the linearized fully compressible formulation confirm the accuracy of the approach. Then, it is subsequently used to characterize the range of existence of instabilities as a function of the Reynolds number, the Mach number, the aspect ratio of the hole, and (for the cavity configuration) the dimensionless volume of the cavity.

Introduction

Plates with orifices are very common elements adopted in numerous industrial applications, like for example silencers, fuel injectors or wind instruments. Under the effect of a harmonic incident acoustic wave, the vortex sheet formed at the lip of the aperture becomes periodically modulated and acts as an amplifier due to a Kelvin-Helmholtz instability, reorganising the jet into an arrangement of vortex rings. The generation of vorticity is an efficient mechanism to dissipate acoustic energy, as a consequence, the use of multiple perforated plates traversed by a mean flow is widely employed as a sound attenuator device for industrial applications, such as gas turbine combustion systems. These systems may suffer from thermoacoustic instabilities because of the potential for unsteady heat release, which can damage the combustion system itself. The flow through a perforated liner with bias flow has been studied experimentally by [START_REF] Heuwinkel | Experimental investigation of the acoustic damping of perforated liners with bias flow[END_REF] while [START_REF] Hughes | The absorption of sound by perforated linings[END_REF]; [START_REF] Eldredge | The absorption of axial acoustic waves by a perforated liner with bias flow[END_REF]; [START_REF] Rupp | The use of perforated damping liners in aero gas turbine combustion systems[END_REF] have conducted both experimental and theoretical investigations. There are, however, situations where the flow through a hole can lead to the opposite effect, namely spontaneous self-oscillations and sound emissions. A particularly favourable situation with respect to sound emission is the flow through two successive holes, as encountered for instance in bird calls and tea kettle whistles [START_REF] Henrywood | The aeroacoustics of a steam kettle[END_REF][START_REF] Longobardi | Studying sound production in the hole-tone configuration using compressible and incompressible global stability analyses[END_REF]. Although less common, the flow through a single hole can also lead to powerful sound emissions. As in other related examples of aeroacoustic resonators, (including for instance the "edge-tone" encountered in the mouthpiece of a recorder or organ pipe for instance), two situations may occur. In the first one, the frequency of the whistling may be directly selected by that of an acoustic resonator located in the vicinity of the hole. This is for instance the case for the so-called "pipe tone" (or pfeifenton), corresponding to a long cylindrical pipe terminated by an aperture of smaller section. In this configuration, which was intensively investigated experimentally by Anderson (1954), the frequency of the whistling directly corresponds to one of the resonance frequencies of the pipe. In the second one, the frequency may be selected by the flow itself regardless of the existence of any acoustic resonator. This situation was noted by [START_REF] Bouasse | Instruments à vent[END_REF], who observed that the flow through a hole in a thick plate separating two large chambers leads to a whistling with a frequency proportional to the thickness of the hole. This observation was rediscovered by [START_REF] Jing | Effect of plate thickness on impedance of perforated plates with bias flow[END_REF]; [START_REF] Su | Measurements and computational fluid dynamics predictions of the acoustic impedance of orifices[END_REF] who, in an effort to improve the performance of perforated plates as sound dampers, reported that, in some circumstances, these devices could lead to self-sustained whistling. In music acoustics one observes the interaction between the two type of mechanisms cf. [START_REF] Coltman | Jet drive mechanisms in edge tones and organ pipes[END_REF]). In the case of the flue instrument, the so-called edge-tone oscillation can coexist with the pipe-tone and under some specific circumstances, as for example during the attack transients, it may be dominant, cf. [START_REF] Castellengo | Acoustical analysis of initial transients in flute like instruments[END_REF]. [START_REF] Verge | Sound production in recorderlike instruments. ii. a simulation model[END_REF] proposed a lumped model for flue instruments where these two feedback loops can coexist and interfere: a hydrodynamic loop responsible for the edge-tone and a cavity loop responsible for the pipe-tone. On the contrary, in the case of the flow past an aperture both mechanisms are associated with the same feedback loop, which is modified by placing a cavity upstream of the perforation. These two situations respectively correspond to the so-called class III and class II categories of aeroacoustic resonators, following the classification of [START_REF] Chanaud | Aerodynamic whistles[END_REF].

Recently, [START_REF] Fabre | The acoustic impedance of a laminar viscous jet through a thin circular aperture[END_REF] used the LNSE approach to investigate the unsteady flow through a circular aperture in a thin plate subjected to harmonic forcing. A novel non-reflecting boundary condition called complex mapping method [START_REF] Sierra | Efficient stability analysis of fluid flows using complex mapping techniques[END_REF] was introduced to overcome the numerical difficulties created by the strong spatial amplification of the fluctuations. The approach allows computing in a rigorous way the impedance of the hole, namely the ratio between unsteady pressure difference across the orifice and unsteady volume flow rate through the orifice, a quantity which can be directly introduced in more elaborate acoustical models. In that study, the authors confirmed that the LNSE can be effectively adopted to predict the impedance even in cases where the spatial evolution of the perturbations is rapidly dominated by nonlinear effects. The same approach was subsequently used by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] for the case of a hole through a thick plate. An important result is that for sufficiently thick holes the impedance can acquire a negative real part in some ranges of forcing frequencies, indicating that energy can be extracted from the flow, thus providing a source for selfoscillations. Investigation of the structural sensitivity also allowed to demonstrate that the hydrodynamic instability of the shear layer separating the jet from the recirculation bubble is the driving motor for the observed phenomenon. This corresponds to the same instability as in the jet of a flue instrument or the shear layer for a grazing flow along a cavity, cf. [START_REF] Dai | Acoustic of a perforated liner with grazing flow: Floquetbloch periodical approach versus impedance continuous approach[END_REF], 2018). In this flow configuration, the sharpness of the aperture corner creates a recirculation bubble that enhances the instability mechanism.

The response of a system to a harmonic forcing is naturally studied via a transfer function: here it corresponds to the concept of impedance, which can also be used to obtain important information regarding the stability properties of a system. First, plotting the impedance in the form of Nyquist diagrams (namely a parametric representation of Z(ω) in the complex plane for real values of ω) provides a direct way to determine the number of unstable modes of the system, as a function of the number of times the Nyquist contour encircles the origin. Secondly, when the system has a complex eigenvalue located close to the real ω-axis, an approximation of the eigenvalue can be obtained from a Taylor expansion of the impedance function around the real axis. Such methods are widely used in several fields like in automatics or in electronics, but remain underemployed in the flow instability community where eigenvalue computation remains the preferred approach. Note, however, that the second idea was recently applied successfully by [START_REF] Sabino | Vortex-induced vibration prediction via an impedance criterion[END_REF] for the problem of vortex-induced vibrations for a spring-mounted cylinder.

The links between impedance and stability properties were explored by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] for the jet flow through a hole. The discussion revealed the existence of two different instability mechanisms leading to sound production: a purely hydrodynamic instability characterised by spontaneous self-oscillations existing in the absence of any incoming acoustic wave, and a conditional instability due to an over-reflection of acoustic waves. Simple criteria formulated in terms of the impedance were given for both kinds of instabilities, allowing to determine their range of existences as function of the hole aspect ratio and the Reynolds number. Among the studies considering a multiply perforated plate, [START_REF] Jing | Effect of plate thickness on impedance of perforated plates with bias flow[END_REF]; [START_REF] Su | Measurements and computational fluid dynamics predictions of the acoustic impedance of orifices[END_REF] measured experimentally the impedances for several configurations with variable hole thickness parameter values, which are in good accordance with the first branch of conditional unstable modes, cf Fabre et al. (2020, Sec 8.). [START_REF] Moussou | An acoustic criterion for the whistling of orifices in pipes[END_REF] studied experimentally a long pipe with a constriction for a number of values of the constriction ratio and the thickness ratio. In this study, they identified both the first and second branch of conditionally unstable modes.

In the approach of [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF], the flow was assumed as locally incompressible, a hypothesis which is expected to be valid for small values of the square of the Helmholtz number (He 2 = ω 2 M 2 ), and which does not directly allow predicting the acoustic field. Nevertheless, they suggested that the locally incompressible solution could be matched to outer solutions incorporating compressibility effects, leading to more elaborate models applicable in situations incorporating for instance acoustic resonators and radiation in an open domain.

The object of the present paper is precisely to show how the impedance computations of [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] based on locally incompressible solution can be used to build a model applicable for a realistic situation involving compressibility. In addition, acoustic pressure fields, obtained from full compressible LNSE computations complement the study. Two generic situations are considered. In the first situation referred to as cavity/open configuration, the domain located upstream of the hole is considered as a closed cavity of finite volume, while the downstream domain is considered as open. We then show that the presence of the upstream resonator can effectively lead to instabilities, as predicted by the conditional instability criterion. The second situation, referred to as open/open configuration, corresponds to the case where the two regions, upstream and downstream f the hole, are considered as open domains of large dimension. We show, in this case, that an instability of purely hydrodynamic type can arise.

The paper is organized as follows. In Sec. 2, the two generic situations are introduced, and the parameters are outlined. In Sec. 3, we show introduce an asymptotically matched or lumped model which allows defining a global impedance for the selected configuration by combining the hole impedance as computed by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] and the impedances of the upstream and downstream domains. We also show that a Taylor expansion of this impedance around the real ω-axis can be used to obtain an instability criterion and an estimation of the eigenvalue of the unstable modes in the fully compressible case. In Sec. 4, we introduce a numerical resolution method for the eigenvalue problem in a fully compressible setup. In Sec. 5 we present results for the cavity/open configuration. We compare both approaches, demonstrating that the asymptotic model is effectively accurate for low Mach numbers. We then provide a parametric study for both problems, thanks to the asymptotic model. Section 6 presents results for the open/open configuration. We particularly investigate the effect of compressibility on the purely hydrodynamic instability mechanism identified by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF], and also consider the acoustic directivity of far-field sound emission.

Problem definition

Fluid parameters

The fluid is considered as a perfect gas with specific constant R g and adiabatic index γ = 1.4. We denote with ρ 0 the reference density and with T 0 the reference temperature (both corresponding to the values in the upstream domain). The fluid is assumed to have constant dynamic viscosity µ and heat conductivity α. The mass flow rate at the inlet of the domain is denoted with ṁ0 , while the mean velocity across the hole is U M = ṁ0 /(ρ 0 πR 2 h ). Based on this velocity scale and the hole diameter D h = 2R h , the Reynolds and Mach numbers of the flow are then defined as:

Re = ρ 0 D h U M µ ≡ 2 ṁ0 πR h µ ; M = U M c 0 with c 0 = γR g T 0 .
(2.1)

The fluid is also characterized by a Prandtl number P r = ρ 0 α/µ which is here assumed to be P r = 0.7. By "open domain" we mean that acoustic waves generated at either side of the hole propagate towards infinity without reflection. Denoting with L h the thickness of the plate, the geometry is thus completely defined by a single dimensionless parameter, the aspect ratio of the hole, defined as:

β = L h 2R h = L h D h (2.2)
In the fully compressible simulations, boundary conditions have to be applied at the boundary of the domain. For simplicity, a half-spherical boundary is considered upstream, and a uniform radial velocity is imposed, as sketched in the figure. Non-reflective boundary conditions used for the compressible computations are introduced and explained in details in section 4.

Cavity/open configuration

The second considered configuration, termed cavity/open configuration, is sketched in fig. 2. This configuration is selected here to study, in the simplest possible setting, the coupling of the hole with a cavity acting as a resonator. The upper domain is considered as a cavity of dimensions L in , R in which acts as a Helmholtz resonator. Therefore, only its volume is relevant, not the exact dimensions L in , R in or the particular geometry. Thus, in addition to the aspect ratio β defined above, a second geometrical parameter enters the problem, namely the dimensionless volume defined as:

V in = L in πR 2 in R 3 h .
The inlet condition is imposed at the leftmost boundary where, for simplicity, a constant velocity profile is enforced, as sketched in the figure.

Matched Asymptotic model

Before considering the resolution of the problem in a fully compressible setting, we detail here a matched asymptotic model which allows us to compute a total impedance characterising the behaviour of linear perturbations of the full system. We first explain how the different regions of the flow domain can be described to obtain the model, and then discuss how the derived total impedance can be used to predict the onset of instabilities.

Matching principle

Under the hypothesis that the Mach number is small and that acoustic wavelengths are much larger than the dimensions of the hole (acoustic compactness hypothesis), it is possible to assume that the flow in the vicinity of the hole is locally incompressible, while compressibility is only relevant in the upstream and downstream domains. This hypothesis is at the origin of the asymptotically matched or lumped model. The ingredients required for matching are the pressure p in (t) just upstream of the hole, the pressure p out (t) just downstream of the hole, and the volume flow rate q(t) across the hole. Working in the frequency domain, all these quantities are expanded as a constant value associated with the base flow, plus a perturbation with harmonic dependency e -iωt , where ω is a (possibly complex) dimensional frequency:

p in (t) = P in + p ′ in e -iωt , p out (t) = P out + p ′ out e -iωt, q(t) = Q 0 + q ′ e -iωt, .
(3.1)

It is important to understand that here p ′ out corresponds to the level of the fluctuating pressure field at distances ∥x∥ considered large in relation to the hole dimension but small compared to the acoustic wavelength, i.e., R h ≪ ∥x∥ ≪ λ with λ = 2πc 0 /ω. Consequently, this pressure level corresponds both to the outer limit for the inner solution, and to the inner limit for the outer solution of the classic matched asymptotic expansion procedure. The same holds for p ′ in which is also used as a matching limit between inner and outer solutions.

Impedance modelling

Inner region : hole impedance

The inner region, located in the vicinity of the hole (delimited by dotted lines in figures 1 and 2), is governed by the incompressible LNSE. A resolution method for this problem was introduced and validated in [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]: for the benefit of the reader, this approach is also briefly summarised in appendix A. The cited method allows ultimately to deduce the (dimensionless) hole impedance Z h (ω) defined as

Z h (ω) = R 2 h ρ 0 U M p ′ in -p ′ out q ′ (3.2)
Here the factor R 2 h /(ρ 0 U M ) is introduced to turn the impedance into a dimensionless one, since the dimensional impedance (p ′ in -p ′ out )/q ′ has physical units kg • s -1 m -4 in the international system, and it is a function of the dimensionless frequency ω = R h ω/U M .

The impedance is ultimately searched as Z h = P • (LN S + iωB) -1 • F, where F represents a forcing of the LNSE by an imposed flow rate, (LN S + iωB) -1 is the linear resolvent of the incompressible LNSE, and P is an operator allowing to extract the overall pressure jump from the linear perturbation. After a convenient discretization, computation of the impedance is thus straightforward and only requires inversion of a single linear problem. It is thus much faster in comparison to eigenvalue computation, which, using shift-and-invert method, typically requires numerous iterative resolutions of such problems. Once Z h (ω) is computed and tabulated [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]) a complete parametric study in terms of Mach number and the cavity volume can be performed as shown below.

Downstream region : radiation impedance

When observed from a large distance, the hole can be seen as a monopolar source, which classically gives rise to spherical diverging waves. This is classically described by a radiation impedance defined as the ratio between pressure p ′ out and flow rate q ′ . This impedance can be obtained by asymptotically matching an acoustically compact inner solution with a monopolar acoustic source, cf [START_REF] Fletcher | The physics of musical instruments[END_REF][START_REF] Pierce | Acoustics: an introduction to its physical principles and applications[END_REF][START_REF] Rossing | [END_REF]). The computation is also reproduced in Fabre et al. (2020, App. A). When expressed in dimensionless variables, the result is a purely real impedance Z rad given by

Z rad = R 2 h ρ 0 U M p ′ out q ′ = M ω 2 2π 3.2.

Upstream region : case of an open domain

In a similar way, in the case of the "open domain" (fig. 1), one can introduce the impedance of the inlet domain Z in , which is defined as

Z in = R 2 h ρ 0 U M p ′ in q ′ = - M ω 2 2π = -Z rad .

Upstream region : case of a cavity

In the case where the upstream domain is considered as a closed cavity (fig. 2), we assume that this cavity acts as a Helmholtz resonator, namely the pressure p ′ = p ′ in , and the density ρ ′ = ρ ′ in are uniform. Then a mass budget leads to

V in c 2 0 dρ ′ in dt = V in dp ′ in dt = ρ 0 q ′
which allows the introduction of the impedance of the cavity Z cav

Z cav = R 2 h ρ 0 U M p ′ in q ′ = i ωM 2 V in = i ωχ , χ = M 2 V in
Note that this expression indicates that the cavity acts as a capacitor for an electrical circuit or as a spring in a mechanical system. Moreover, its characteristics only depend upon the quantity χ = M 2 V in which combines the Mach number and the dimensionless volume of the cavity. Such a model could be complemented with the addition of two other terms that have been neglected. In particular, one could include on the left-hand side of the mass balance the deviation from isentropic pressure fluctuations due to, for instance, the effects of the thermal boundary layer and on the right-hand side, the convective term involving density fluctuations. These terms have been neglected based on the fact that velocity and temperature gradients within the cavity are small as long as the ratio between the height of the cavity and the radius of the hole is large. In our study, this corresponds to L in = R in = (V in /π) 1/3 ≫ R h , which holds for the cavities analysed in this study.

Summary : total impedance of the problem

Regrouping all the regions, we are now able to obtain a single constitutive equation for the total impedance of the system, denoted either Z a or Z b for the two investigated configurations, which allows us to determine the eigenfrequencies of the complete problem:

Z a (ω) = Z h (ω) + M ω 2 π = 0 (3.3) (b): for the cavity/open configuration, Z h = -Z cav -Z rad , or equivalently : Z b (ω) = Z h (ω) + M ω 2 2π + i M 2 V in ω = 0 (3.4)
We emphasize that the total impedance defined here is designed mainly to be used to detect eigenvalues, hence only the condition Z(ω) = 0 is significant. The complex zeros ω = ω R + ω I of Z then correspond to the eigenmodes of the system, and the system is therefore unstable if there exists such a zero with ω I > 0. When not zero, there is no direct physical interpretation to the value Z(ω) associated with a given ω. Schematically, 1/Z it can be conceived of as a measure of the response of the system to an imposed forcing, so that Z = 0 means that the response is infinite, or in other words that a solution without forcing is possible. For instance, Fabre et al. (2020, App. A) considered the case where the forcing corresponds to a spherically converging wave coming from downstream; in this case the reflection coefficient is effectively proportional to Z -1 (see Eq. A12 in this reference). Other kinds of forcing could be considered, leading to the same conclusion. In the present paper we remain to an intuitive interpretation of Z -1 and do not elaborate on the link between the impedance and any specific forcing.

Predicting instability from a Taylor expansion of the impedance

As stated in the introduction, knowledge of the impedance function Z(ω) along the real ω-axis allows obtaining important information regarding instability properties of the system in two ways. First, Cauchy's argument principle can be used as a graphical method to determine whether or not an instability exists. This argument is developed in appendix C. Second, eigenvalues located close to the real axis may be expected to be accurately predicted from a Taylor expansion of the impedance around the real axis. This argument is presented here.

Asymptotic prediction of eigenvalue for the cavity/open configuration

Following an idea previously used in Ferreira [START_REF] Sabino | Vortex-induced vibration prediction via an impedance criterion[END_REF] for the problem of vortex-induced vibrations of a spring-mounted cylinder, we assume here that the impedance of the full system is mostly reactive. In the present case, this means that the impedance is dominated by its imaginary part, while the real parts (i.e., Re(Z h ) and the radiation impedance) correspond to lower order terms. Such a hypothesis, together with the fact that the flow is acoustically compact within the region of the hole, allow the use of an asymptotic expansion truncated at first order to determine the zeros of the impedance. We first elaborate this idea for the cavity/open configuration. The hypotheses are as follows:

(

i) |Re(Z h )| ≪ |Im(Z h )|, i.e. |Re(Z h )| ∼ ε |Im(Z h )|, (ii) M ω 2 2π ≪ |Im(Z h )|, i.e. M ω 2 2π ∼ ε |Im(Z h )|.
where the real parameter 0 < ε ≪ 1. Note that hypothesis (i) is not justified for every value of ω 0 since from the results of [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] the real and imaginary parts of Z h are generally of comparable order of magnitudes. However, this hypothesis can be expected to be valid in the vicinity of the threshold of the instability. Hypothesis (ii) is needed for the acoustic compactness and therefore directly satisfied.

Consider the frequency expansion

ω = ω 0 + εω 1 , ω 0 ∈ R, ω 1 ∈ C, ε ∈ R, (3.5)
and let substitute ω in eq. (3.4) by eq. (3.5) and by performing a Taylor expansion in terms of the assumed small quantities leads to

Z b (ω) = i Im Z h (ω 0 ) + 1 M 2 V in ω 0 +ε Re Z h (ω 0 ) + M ω 2 0 2π + ∂Z h ∂ω ω=ω0 - i M 2 V in ω 2 0 ω 1 +O ε 2 (3.6)
where O ε 2 denotes higher order terms as function of the assumed small parameter. The condition Z b = 0 then leads to the following results:

(i) The zeroth-order terms lead to the condition

-ω 0 Im Z h (ω 0 ) = 1 M 2 V in = 1 χ (3.7) (ii)
The first-order term leads to

Im(ω 1 ) = Re Z h (ω 0 ) + M ω 2 0 2π ∂Im Z h (ω) ∂ω R ω=ω0 - 1 χω 2 0 ∂Re(Z h (ω) ∂ω R ω=ω0 2 + ∂Im Z h (ω) ∂ω R ω=ω0 - 1 χω 2 0 2 , Re(ω 1 ) = -Re Z h (ω 0 ) + M ω 2 0 2π ∂Re(Z h ) ∂ω R ω=ω0 ∂Re Z h (ω) ∂ω R ω=ω0 2 + ∂Im Z h (ω) ∂ω R ω=ω0 - 1 χω 2 0 2 (3.8)
The imaginary part of the first-order correction directly provides a criterion of stability.

Provided that the imaginary part of ∂Z h /∂ω is negative (a condition which is found to hold in all cases where the starting hypotheses are verified), then it is possible to conclude that an instability is possible as soon as

Re Z h (ω 0 ) < - M ω 2 0 2π .
We recognize here an improved version of the conditional instability criterion of [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]. Physically, this condition means that the energy extracted from the baseflow

-Re Z h (ω 0 ) |q ′ | 2 /2 must be larger than the energy radiated Z rad |q ′ | 2 /2 Fabre et al. (2020) have documented the function Z h (ω)
for real values of ω over a wide range of parameters. Once the hole impedance function Z h (ω) is determined, these results can be used to solve the coupled conditions eq. (3.7),eq. (3.8) and ultimately to obtain an instability criterion and an estimate for the growth rate. Figure 3 explains graphically these conditions. The resolution can be done in two ways. Via a direct method, that is, given the parameters M and V in , one first solves for eq. (3.7), which is an implicit equation in ω 0 as a function of the parameter χ (as sketched in fig. 3a). Then one may deduce Im(ω 1 ) which is an explicit function of ω 0 and M (as sketched in fig. 3b) and it ultimately provides a criterion of instability.

An alternative is to follow an inverse method. Given M , we first consider Im(ω 1 ) as a function of ω 0 and deduce the ranges of ω 0 where this function is positive (as indicated in blue on fig. 3b). Once these unstable ranges are known, we deduce the corresponding ranges for 1/(M 2 V in ) by using eq. (3.7) (as indicated by blue ranges in fig. 3a). The approach thus indicates the ranges of V in where, for the given M , the jet is unstable.

The great advantage of this inverse method is that the equation eq. (3.7) is explicit when solving for V in in terms of ω 0 . The inverse method is a very efficient way to obtain an estimation of the eigenvalue of the full problem ω = ω(Re, β, M, V in ) provided one disposes of a tabulation of the function Z h (ω; Re, β) for real values of ω. It must be emphasised that the number of parameters has been reduced from four to only two, as V in and M only occur through the modelled impedance of the upstream and downstream domains. However, the reduction relies on a series of strong hypotheses: first M ≪ 1 and |ω| ≪ 1/M for the matched asymptotic model to hold, and second the assumptions used to treat Re(Z h ) as a correction. The validity of the approach, therefore, have to be assessed by comparing the results with those obtained using a fully compressible model in order to clarify the range of validity of the used approximations, as detailed in section 5.1.

Asymptotic prediction of eigenvalue for the open/open configuration

Let us now follow a similar route to achieve an estimation of the eigenvalue ω for the open/open configuration. In this case, the zeroth-order and first-order corrections simplifies to

-ω 0 Im(Z h )(ω 0 ) = 0 (3.9) Im(ω 1 ) = Re Z h (ω 0 ) + M ω 2 0 π ∂Im(Z h ) ∂ω R ω=ω0 | ∂Z h ∂ω ω=ω0 | 2 , Re(ω 1 ) = -Re Z h (ω 0 ) + M ω 2 0 π ∂Re(Z h ) ∂ω R ω=ω0 | ∂Z h ∂ω ω=ω0 | 2 (3.10)
where the non-zero ω 0 are the zeros of the imaginary part of the hole impedance function Im(Z h ) and the growth rate is estimated by eq. (3.10).

Note that this expression is identical to the one obtained for the cavity/open configuration when V in → ∞, except for the radiation term, which is twice the value in the previous case. This accounts for the fact that radiation occurs on both sides, so that total radiation losses are twice larger.

Full compressible formulation

After detailing the matched asymptotic model, we now introduce in this section a numerical method to resolve directly the eigenvalue problem in a fully compressible setting.

Compressible Navier-Stokes equations

Let's consider a compressible fluid motion of a perfect gas described in primitive variables by q = [ρ, u, T, p]

T , where the velocity vector field is u = (u, v, w), pressure p, temperature T and fluid density ρ. Dimensional primitive variables have been made dimensionless, as follows:

x = x D h , t = tU M D h , ρ = ρ ρ ref , u = ũ U M , T = T T ref , p = p -p ref ρ ref U 2 M (4.1)
where dimensional values are designated by an upper tilde •, and reference values are indicated with the • ref . Dynamics is governed by the compressible Navier-Stokes equations, which are here written in terms of primitive dimensionless variables in the compact vector notation:

M ∂q ∂t = N S q = L q + N(q) + C = 0 (4.2)
where C = [0, 0, 0, 1] T , the mass matrix M and the linear operator L are defined as

M =     1 0 0 0 0 ρI 0 0 0 0 ρ 0 0 0 0 0     , L =     0 0 0 0 0 -∇ • τ (•) 0 ∇ 0 0 -γ P rRe ∆ 0 0 0 0 γM 2     (4.3)
while the nonlinear operator is written as

N(q) =      u • ∇ρ + ρ∇ • u ρu • ∇u (γ -1) ρT ∇ • u -γM 2 τ (u) : D(u) + ρu • ∇T -ρT      , (4.4)

Compressible Navier Stokes -Base flow equations

The stability of a steady-state solution to infinitesimal perturbations can be analysed using the classical approach based on linearization of the governing equations: the total flow field is expanded into the sum of a steady-state term plus an infinitesimally small unsteady harmonic perturbation as follows:

q(t) = q 0 + ε qe -iωt + c.c. (4.5)
where ε ≪ 1. Inserting eq. (4.5) in the governing equations eq. ( 4.2) and neglecting quadratic terms leads to two problems, one for the base flow and one for the perturbation.

In particular, at leading order, only steady terms are kept, which leads to the steady Navier-Stokes equations:

N S q 0 = L q 0 + N(q 0 ) + C = 0 (4.6)
complemented with appropriate boundary conditions. No-slip adiabatic boundary conditions are used at the walls (eq. (4.7c)). At the axis of revolution, the radial component v 0 is set equal to zero, and the radial derivative of the remaining terms is null eq. (4.7d)). At the outlet we set stress-free and isothermal boundary conditions (eq. (4.7b)); in this way the pressure at the outlet is equal to the thermodynamic pressure, i.e., p 0 = 1. Finally, at the inlet boundary Γ in = Γ in,0 ∪ Γ in,1 , a constant mass flow is enforced on the Γ in,0 boundary, slip condition, constant density and zero thermal flux are imposed on the Γ in,1 (eq. (4.7a)). Summarising,

ρ 0 | Γin = 1, Γin,0 ρ 0 u 0 • ndS = π 4 , u 0 • n | Γin,1 = 0, ∇T 0 • n | Γin = 0 (4.7a) T 0 | Γout = 1, -p 0 I + τ (u 0 ) • n Γout = 0 (4.7b) u 0 | Γw = (0, 0, 0) T , ∇T 0 • n | Γw = 0 (4.7c) v 0 | Γw = 0, ∂u 0 ∂r = ∂w 0 ∂r = ∂ρ 0 ∂r = ∂T 0 ∂r =
∂p 0 ∂r = 0 on Γ axis (4.7d )

Linearized compressible Navier-Stokes equations -Homogeneous problem

The linearized compressible Navier-Stokes equations govern the evolution of the perturbation q:

-iωMq = LN S 0 q = L + DN q0 q, (4.8)

where DN q0 is the jacobian matrix of the nonlinear operator evaluated at the steadystate q 0 . With the purpose of modelling a large container upstream of the hole, for the open/open case we have designed a computational domain, fig. 4

ρ| Γin = 0, -pI + τ (û) • n Γout = 0, ∇ T • n | Γin = 0 (4.9a) ρ| Γout = 0, -pI + τ (û) • n Γout = 0, ∇ T • n | Γout = 0 (4.9b) û| Γw = (0, 0, 0) T , ∇ T • n | Γw = 0 (4.9c) v| Γw = 0, ∂ û ∂r = ∂ ŵ ∂r = ∂ ρ ∂r = ∂ T ∂r = ∂ p ∂r = 0 on Γ axis (4.9d )
In particular, in the far-field (inlet and outlet) we impose null density variations, stressfree boundary condition and vanishing thermal flux (eqs. (4.9a) and (4.9b)); doing so the mass flux, ρ 0 û • n, is allowed to vary. A no-slip adiabatic boundary condition is used at the walls (eq. (4.9c)), while at the axis the radial component of the velocity v is set to zero, together with null radial derivative of the remaining terms (eq. (4.9d)).

For the purpose of modelling a closed cavity that acts as an acoustic resonator, we have a computational domain, which is sketched in fig. 4 (b), where the complex mapping layer is only present in the region placed downstream of the hole. The set of boundary conditions are as follows

ρ| Γin = 0, û| Γin = (0, 0, 0) T , ∇ T • n | Γin = 0 (4.10a) ρ| Γout = 0, -pI + τ (û) • n Γout = 0, ∇ T • n | Γout = 0 (4.10b) û| Γw = (0, 0, 0) T , ∇ T • n | Γw = 0 (4.10c) v| Γw = 0, ∂ û ∂r = ∂ ŵ ∂r = ∂ ρ ∂r = ∂ T ∂r = ∂ p ∂r = 0 on Γ axis . (4.10d )
i.e. null density and velocity variations (eq. (4.10a)) at the inlet, stress-free condition, null density variation and vanishing thermal flux (eq. (4.10b)) at the outlet; no-slip and adiabatic walls (eq. (4.10c)); null radial velocity component and null radial derivative of the remaining terms (eq. (4.10d)).

Numerical implementation

Following a usual route in global stability analysis, the nonlinear problem eq. (4.6) for the base flow is solved using a Newton iteration method and the eigenvalue problem eq. (4.8) is solved using a shift-invert Arnoldi method. Spatial discretization is done using a finite-element method, using P2-elements for velocity components u x , u r and P1-elements for thermodynamic variables P, ρ, T . Mesh generation and assembly of matrix operators is performed using the FreeFem++ software [START_REF] Hecht | New development in freefem++[END_REF]. Resolution is achieved using PETSc/SLEPc libraries, which are directly implemented in FreeFem++. Monitoring of computation, loop over the parameters and post-processing are handled in Matlab thanks to the StabFem suite [START_REF] Fabre | A practical review on linear and nonlinear global approaches to flow instabilities[END_REF]. Note that during the process, mesh adaptation is used in a way similar as described in [START_REF] Fabre | A practical review on linear and nonlinear global approaches to flow instabilities[END_REF], to ensure that the resolution is sufficient to ensure grid-independence when computing the baseflow and the eigenmodes. Examples of codes reproducing sample results are shared on the website of the StabFem project.

Results -Cavity/open configuration

Validation of the asymptotic model -Comparison with compressible LNSE

In section 3.3.1, we introduced an asymptotic method which is able to predict the eigenvalues ω = ω(Re, β, M, V in ) from a simple tabulation/computation of the function Z h (ω; Re, β) for real values of ω, hence reducing the number of parameters from four to two only. Before conducting a full parametric study of the instability with the proposed matched asymptotic method, we have to assess its validity by comparing the predictions with resolution of the full eigenvalue problem. This is done in fig. 5 which compares the amplification rates (left column) and frequencies (right column) obtained with the two approaches for values of M , Re and V in spanning a large range of parameter values, considering a hole with aspect ratio β = 0.3.

Consider, first, the predictions of the asymptotic model represented by coloured lines in the figures. Thanks to the inverse method explained in section 3.3.1, the asymptotic prediction allows us to plot ω as a continuous function of V in . We use solid lines for the segments of the curves corresponding to unstable modes and dotted lines for those corresponding to stable modes. As identified by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF], for β = 0.3, several modes of conditional instability, termed C 1 , C 2 , etc. . . are expected to arise as the Reynolds number is increased. The corresponding frequencies are quantized, and an argument to explain this quantification was proposed in terms of the dynamics of the shear layer. An alternative argument, in terms of a forward shear wave and a backward acoustic wave, is proposed in appendix D.

From the results of [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF], for β = 0.3 (see also fig. 18 reproduced in appendix D), the first mode C 1 arises just below Re = 800 and the second mode C 2 arises for Re ≈ 1500. This is in good agreement with the observed results of the asymptotic model, which effectively predicts two ranges of instability for Re = 1600 and Re = 2000, at least, for the smallest considered values of M . The figure also shows that increasing M results in a shifting of the instability ranges towards smaller values of V in .

Consider now, the eigenvalue calculations, represented by circles in the figures. Results have been computed for a limited set of values of V in where unstable modes were expected. Recall that in the eigenvalue study, V in is linked to the size of the numerical domain, so that the whole process (mesh generation, base flow computation, resolution of eigenvalue problem) has to be restarted for each new value of V in . An excellent matching between the two estimates may be appreciated even for large growth rates, the relative error being less than 3% in most cases. Comparison seems poorer, at first sight, for the case of Re = 800 reported in fig. 5 (a) but it must be reminded that the case is very close to the threshold and amplification rates are very small, so that the absolute error is actually of comparable order to the other cases. Not that an excellent agreement is still found in cases where the amplification rate is not small, a range where the impedance criterion should be a priori slightly less reliable due to its perturbative nature. The agreement also remains excellent when the Mach number is raised to M = 2 • 10 -2 . Note that for eigenvalue computations it has been only considered configurations with V in > 10 2 , since for smaller values the cavity becomes very small and the modelling as a Helmholtz resonator becomes questionable. This is why we did not attempt to draw any comparisons for M > 2 • 10 -2 , with the exception of a case with M = 4 • 10 -2 represented in plot (e).

Structure of some eigenmodes

Let us now illustrate the structure of a few eigenmodes computed with the full compressible LNSE. Figure 6 displays the eigenmode computed for M = 5 • 10 -3 and V in = 10 4 for Re = 1200. This mode is correctly predicted by the asymptotic model, and recognized to correspond to the branch C 1 of conditional instability modes, as defined by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]. As observed, the pressure level inside the cavity is uniform, confirming that the cavity effectively acts as a Helmholtz resonator for this mode so that the modelling hypotheses are correctly verified. Downstream of the hole, the mode is characterized by an alternance of structures of opposite sign, localized along the shear layer. This structure is characteristic of regions associated with negative real part of impedance, as identified by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]. Note that very far away in the downstream domain, the structure is expected to match with a spherical diverging wave of the dimensionless wavelength λ = 2π/(M ω R ) . Here λ is of order 70, so this structure is not visible on the figure. A characterization of the far-field acoustic radiation is described in section section 6.3, see fig. 13.

In addition to eigenmodes of the kind presented in Figure 6 which are well predicted by our asymptotic approach, one typically observes the existence of other families of eigenmodes with a more complex structure. Figure 7 (b) displays a family of such modes, computed for the set of parameters M = 2 • 10 -1 , V in = 10 4 for Re = 1200. One clearly observes that the pressure inside the cavity is no longer uniform, but characterised by nodal lines in the radial and axial distributions. These modes are recognised as cavity modes. They arise as soon as the acoustic compactness hypothesis fails, i.e., when the acoustic wavelength is smaller than the characteristic length (L in = (V in /π) 1/3 ) of the cavity.

Parametric study

In our previous work, the ranges of parameters corresponding to a conditional instability (requiring the presence of a correctly tuned resonator) were mapped in the Re -β plane : see Fabre et al. (2020, Fig. 13), also reproduced in appendix D (see fig. 18). We are now able to build upon these results a parametric study of the situation where the resonator corresponds to the upstream cavity, as function of the four parameters (Re, β, M, V in ). Figures 8 and9 displays the dependence of the neutral curves on the 9 (a,d). These regions of instability grow with the Reynolds number, and they shrink with V in . Contrary to instabilities C 1 and C 4 , for a given value of V in instabilities C 2 and C 3 may exist for every M , for such a reason these instabilities may be conceived as a degenerate situation of pure hydrodynamic instabilities H 2 and H 3 , which are discussed in section 6. Finally, the dependence of this type of instability on the acoustic resonator is better appreciated if we consider the effect of V in and M together with the χ = V in M 2 parameter. This allows us to display neutral curves of stability in the (χ, Re) plane, Figure 11 (a) displays the range of existence of instabilities as function of Re and β for three values of M . Two different modes, corresponding to the purely hydrodynamic modes H 2 and H 3 as identified by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF], are documented. Note that higher instability modes, called H i with i = 4, 5, 6 . . . , exist in the studied interval of the parameter β ; however, they arise at larger values of Re and are not considered. The curves displayed in the figure for the smallest value of M , namely 10 -2 , are very close to the predictions obtained by Fabre et al. (2020, Figure 15), which is represented by dots in fig. 11 . One can see that compressibility has almost no effect on the instability threshold of mode H 2 . On the other hand, it has a destabilising effect on mode H 3 . Figure 11 (b) investigates the effect of compressibility on the oscillation frequency, here represented as a Strouhal number St β = ω R β 2π . One can see that compressibility decreases the frequency for the shortest holes and has almost no effect on frequency for longer ones. This behaviour is associated with the significant modification of the threshold of instability (Re) for short holes (left asymptote of either H 2 or H 3 in fig. 11), which does not occur for holes with larger β value. A substantial variation in the critical Reynolds induces a modification of the vena contracta coefficient α vena , see e.g., fig. 16, which in turn may be linked to the frequency of the instability, see the discussion in appendix D.

Effect of Mach number -Sensitivity analysis

In order to explain these observed trends, we consider the sensitivity dω dM of the complex frequency ω with respect to M . This quantity may be split into two terms, corresponding respectively to the sensitivity to base flow modifications and to the sensitivity to a Mach number variation of the linearized equations:

dω dM M = ∂ω ∂M q0 + ∂ω ∂q 0 M ∂q 0 ∂M (6.1)
We have employed two techniques, a continuous adjoint technique described in [START_REF] Meliga | Effect of compressibility on the global stability of axisymmetric wake flows[END_REF] and a forward evaluation of the sensitivity. Provided the Mach number is close to the incompressible limit, the continuous adjoint technique provides less accurate results, for such a reason we have decided to perform a forward evaluation. The first term ∂ω ∂M is evaluated using a first order finite difference, which requires the resolution of two eigenvalue problems eq. (4.8) with the steady-state frozen at the Mach number M , and with a Mach number perturbation of a small magnitude ∆M ≪ M in the linearized Navier-Stokes operator of the eigenvalue problem. The second term is also evaluated by finite difference issued of the resolution of the eigenvalue problem at fixed Mach number but with two different steady states computed at M and M + ∆M . Figure 12 displays the sensitivity computed in such a way, for a value of Re corresponding to the thresholds of the instability for modes H 2 (with β = 0.8) and H 3 (with β = 2). The figure also displays the value of ∂ω 1 /∂M obtained from the asymptotic model eq. (3.8), which is a constant. Figure 12 (a) reports a linear variation of the two terms of the sensitivity with respect to Mach number. A Mach number variation in the base flow has a stabilizing effect, whereas the instability is triggered by small variations of the Mach number of the linearized operator. The most dominant term for this kind of acoustically compact solutions seems to be the base flow effect, which has a small stabilizing effect. In particular, it explains the almost insignificant variation of the H 2 neutral curves in fig. 11. Concerning the frequency, both terms have an almost opposite effect, which implies an almost null variation of the instability frequency for the H 2 mode. The impedance criterion in this case predicts a stabilizing effect with small frequency increase, which holds relatively well for M < 0.02. On the other hand, the variation with respect to the Mach number of the growth rate and the frequency of instability of a configuration that is no longer acoustically compact for Mach numbers of the order of M ∼ 10 -2 such as the H 3 mode at threshold for β = 0.8 greatly differs with the estimations made with the impedance criterion. Figure 12 (b) reports a similar stabilizing effect of the base flow to the one of the H 2 mode. However, in this case variations of the Mach number in the linearized operator greatly destabilize the steady-state, which causes the large variations in the neutral curves displayed in fig. 11. In terms of frequency variations, it is possible to observe much larger excursions, which are negative and constant for M < 0.07 and increase linearly for M > 0.07. So we may conclude that the impedance criterion holds relatively well for large β and instability modes with low frequency, which are in turn the most acoustically compact, but it fails to predict accurate trends even for low Mach numbers for modes with higher frequencies and small length to thickness ratios.

Directivity of acoustic emission

Finally, we address the influence of parameters (β, M ) on the directivity pattern of instabilities of type H 2 and H 3 . For that purpose, we evaluate the set of neutral eigenmodes for each pair (β, M ). Notice that the amplitude of the eigenmodes are arbitrary, the pressure levels displayed in fig. 14 have been normalized with respect to the monopole radiation (based on the oscillating volume flux through the perforation). Three values of the Mach number M = {10 -2 , 2 • 10 -2 , 5 • 10 -2 } and two of the dimensionless parameter β = {1, 2} are selected for this study. The configuration β = 1 corresponds to a configuration less acoustically compact than β = 2 and it seems a priori more likely that the radiation differs from the single monopole pattern. Figure 13 displays the

Conclusion

The objective of the present paper was to investigate how the instability potential of a single jet passing through a hole in a thick plate, recently identified by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] using linearized Navier Stokes equations (LNSE) in a strictly incompressible setting, manifests in a more realistic configuration involving compressibility. For this sake, we considered two generic situations. In the first situation, the upstream domain is a closed cavity and the downstream domain is an open space. This situation was chosen to check the conditional instability mechanism, requiring the existence of a conveniently tuned resonator. In the second situation, the two regions, upstream and downstream of the hole, are considered as open. This situation was chosen to check the purely hydrodynamical instability which is expected to exist even in the absence of a resonator.

The two cases, has been analyzed with an asymptotic method, which provides an instability criterion and an estimate for the amplification rate. The method consists, in a first step, in writing an impedance of the global system incorporating the hole impedance as computed by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF]) and modelling the impedance of the upstream and downstream regions, and in a second step in performing a Taylor expansion of this impedance around the real ω-axis to identify its zeros. The great advantage of the method is that the Mach number M and the cavity volume V in appear as parameters in the model, so that a parametric study of the problem can be done entirely in terms of the Reynolds number and aspect ratio, therefore reducing the number of computational parameters from four to two. The potential to accurately predict the instability properties via the asymptotic model is put into test in section 5.3 and section 6 for the conditional and pure hydrodynamic cases, respectively. A cross-comparison with the results carried out with the compressible Navier-Stokes equations shows a good matching between the two approaches. The impedance criterion has been employed to identify the regions of existence in the (Re, χ) plane of a series of C i , i = 1, 2, 3, 4, . . . modes. In addition to these acoustically compact modes, at larger Reynolds numbers, there exist unstable modes associated with higher order modes of the cavity connected upstream (see fig. 7 for an example of that phenomenon). The use of the impedance criterion for the characterisation of the compressibility effect in the pure hydrodynamical case is less accurate. We have observed that the estimations of the growth rate are relatively acceptable for the H 2 mode, but they are faulty for the H 3 mode, in particular for small length to diameter ratios. The inadequacy of the criteria to characterise this case can be attributed to the lack in the asymptotic model of the effect of the backward-travelling acoustic wave. These results suggest that a better modelling of the hydrodynamic-acoustic interaction is required to gain in accuracy. Finally, in section 6 we have examined the influence of the Mach number on the directivity pattern of the family of the pure hydrodynamical modes.

Declaration of interests. The authors report no conflict of interest.

Appendix B. Properties of the compressible steady-state

As illustrated in figure fig. 15, the base flow is characterized by a recirculation region originating from the upper corner of the hole. The pressure jump due to this recirculation can be represented by the so-called discharge coefficient (also called vena contracta coefficient

α vena = R 2 h R 2 J
where R J is the effective radius of the jet. This function has been tabulated by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] as function of Re and β in the incompressible case. When taking into account the compressibility effects, the discharge coefficient [START_REF] Bragg | Effect of compressibility on the discharge coefficient of orifices and convergent nozzles[END_REF]) can be written in term of the dimensionless variables introduced in section 4.1 as

α vena = ṁ0 πR 2 h 1+γM 2 Pin √ Tin 1 M 2 γ-1 ( 1+γM 2 Pin 1+γM 2 Pout ) -2/γ -( 1+γM 2 Pin 1+γM 2 Pout ) -(γ+1)/γ
, which in the low Mach number limit can be approximated as

α vena = ṁ0 πR 2 h ρ in (2(P in -P out ) -3M 2 (P in -P out ) 2
and it coincides with the one employed by [START_REF] Fabre | Acoustic impedance and hydrodynamic instability of the flow through a circular aperture in a thick plate[END_REF] at the incompressible limit. Note that at large Reynolds numbers other simpler estimates exist, for instance see the discussion by [START_REF] Gilbarg | Jets and cavities[END_REF] for the compressible Borda tube, which has been revisited by Durrieu et al. (2001) and compared against experimental evidence. Figure 16 (a) displays the effect of the Mach number on the discharge coefficient. It shows a good quantitative agreement with the theoretical estimation, and it weakly increases with the Mach number. Compressibility effects accelerate the bulk flow within . In addition, it has been observed that the shear layer thickness remains unchanged with a weak Mach number increase (M < 0.2). In section 6 it is shown that an increase of the Mach number in the steady-state solution has a stabilising effect, which can be attributed to an attenuated recirculation region.

Appendix C. Nyquist curve -Cauchy's argument principle

Let us review the use of the Nyquist criterion together with the drawing of Nyquist curves. In the absence of poles in the real axis, the Nyquist plot is drawn along the real axis. However, in the presence of a pole of the impedance in the real axis, one must take a contour that does not encircle the pole. In particular, the augmented system impedance possess a pole at ω = 0, therefore a complex contour that does not encircle zero is employed as the one depicted in fig. 17 (a). In the evaluation of the impedance, let us consider here without loss of generality the augmented system impedance Z (a) (either Z (a) = Z a or Z (a) = Z b ) along the contour Γ , i.e., Z (a) (Γ ), which provides a direct evaluation of the stability of the system. Provided that the contour of integration does not encircle any pole of the system, which is satisfied by construction, the number of times that the curve Z (a) (Γ ) encircles the origin in the counter-clockwise direction determines the number of zeros in the area surrounded by the contour Γ . In the condition that the contour of integration Γ encloses the unstable complex plane, then any encirclement of the origin implies that the system is unstable. To illustrate this, let consider the Nyquist curve represented in fig. 17 (b), where the curve Z (a) (Γ ) is oriented counter-clockwise, and it encircles twice the origin: this implies that the system has two unstable zeros. A more careful evaluation reveals that this corresponds to a pair of conjugated complex zeros.

Additionally, impedance values for real ω can also provide an estimation of the complex zeros whenever Im(ω) is of small magnitude. Here, let us detail the procedure followed in section 3.3.1. We consider the case where the Nyquist curve is found near the origin. The first scenario corresponds to a complex frequency ω = ω R + iω I = ω 0 + εω 1 , where

= Z (a) (ω) = Z (a) (ω 0 ) + dZ (a) dω ω=ω0 (ω -ω 0 ) + R 1 (ω), ( C 1)
where R 1 (ω) is the remainder of the Taylor expansion. The remainder of the expansion can be shown to be bounded, for instance, using Cauchy's integral formula and the maximum principle (see, for instance, [START_REF] Rudin | Complex and real analysis[END_REF]) yields

|R 1 (ω)| ⩽ M r |ω -ω 0 | 2 r(r -|ω -ω 0 |) ⩽ M r η 2 1 -η , with M r ≡ max |ω-ω0|=r |Z (a) (ω)|
for |ω -ω 0 | < r and |ω-ω0| r ⩽ η < 1, where we have assumed that the impedance function is holomorphic in a closed disk of radius r of the complex plane. Therefore, the function R1(ω) ω-ω0 is also analytic within the disk. In such a way we can approximate the value of zero ω as

ω-ω 0 = Z (a) (ω 0 ) dZ (a) dω ω=ω0 + R 1 (ω) ω -ω 0 ≈ Z (a) (ω 0 ) dZ (a) dω ω=ω0 = Z (a) R (ω 0 ) + iZ (a) I (ω 0 ) ∂Z (a) R (ω) ∂ω R ω=ω0 + i ∂Z (a) I (ω) ∂ω R ω=ω0 , ( C 2
) where the error of the approximation is

ϵ ⩽ Z(ω 0 ) dZ dω ω=ω0 η ′ 1 -η ′ with η ′ = |R1(ω)|
|ω-ω0| /| dZ dω ω=ω0 |, because the radius of convergence of the rational complex function 1 dZ dω ω=ω 0 +z is equal to | dZ dω ω=ω0 |. Thus, the approximation eq. (C 2) converges uniformly far from the critical points of the impedance. Furthermore, it provides good estimates whenever |Z(ω 0 )| ∼ ε ≪ 1 and |ω -ω 0 | ∼ ε ≪ 1, which is the motivation to the condition (i) and the expansion in frequency in section 3.3.1. Finally, note that eq. (C 2) could be used as a step in a Newton iteration, whenever the initial guess ω 0 is far from the zero of the impedance.

Multiplication by the complex conjugate of the denominator in eq. (C 2) leads to

εω 1 ≈ Z (a) R (ω 0 ) + iZ (a) I (ω 0 ) ∂Z (a) R (ω) ∂ω R ω=ω0 -i ∂Z (a) I (ω) ∂ω R ω=ω0 ∂Z (a) R (ω) ∂ω R ω=ω0 2 + ∂Z (a) I (ω) ∂ω R ω=ω0 2 . ( C 3)
Finally, one could split eq. (C 3) in real and imaginary parts, this yields the following expression for ω

R εRe(ω 1 ) = ω R -ω 0 ≈ Z (a) R (ω 0 ) ∂Z (a) R (ω) ∂ω R ω=ω0 + Z (a) I (ω 0 ) ∂Z (a) I (ω) ∂ω R ω=ω0 ∂Z (a) R (ω) ∂ω R ω=ω0 2 + ∂Z (a) I (ω) ∂ω R ω=ω0 2 . (C 4) and for ω I εIm(ω 1 ) = ω I ≈ Z (a) I (ω 0 ) ∂Z (a) R (ω) ∂ω R ω=ω0 -Z (a) R (ω 0 ) ∂Z (a) I (ω) ∂ω R ω=ω0 ∂Z (a) R (ω) ∂ω R ω=ω0 2 + ∂Z (a) I (ω) ∂ω R ω=ω0 2 . ( C 5)

Appendix D. Frequency selection argument

This appendix provides an argument explaining the quantification of the eigenvalues observed for the cavity/open configuration.

The frequency of the sound generated is selected by considering the two elements that compose the feedback loop of a Class III aerodynamic whistle: the hydrodynamicacoustic wave interaction and the acoustic resonator. In this kind of mechanism, we can distinguish two feedback loops, a first loop composed of the interaction of a hydrodynamic instability with an acoustic wave and a second one which accounts for the interaction of the first feedback loop with the acoustic resonator. A hydrodynamic-acoustic wave interaction develops whenever the shear layer of the jet is unstable and the jet acts as a source of energy, which occurs when the resistance of the hole is negative. The shear layer instability is triggered at the leftmost corner of the hole, disturbances grow along the hole, however in the case when the shear layer is not sufficiently unstable (so selfsustained oscillations arise by pure hydrodynamical arguments) it requires an acoustic wave to close the loop, which is an instantaneous process in the low Mach number limit. For acoustically compact source regions, the frequency selection of this mechanism is dominated by the hydrodynamic instability, because the travel time of the acoustic wave is of lower order of magnitude. Secondly, the cavity acts a resonator, selecting a set of discrete frequencies among those associated with a negative resistance.

In the present configuration, there exist four branches of instability, each of them denoted as C n for n = 1, 2, 3, 4, which are characterised by a nearly constant Strouhal number: + ∆St β , where ∆St β ≈ [0.6, 0.7]. This frequency shift may be estimated if one realizes that the global instability is the result of the constructive interaction between two travelling waves: a downstream-travelling wave, which is excited about the hole lip and propagates around the jet core position, and an upstream travelling wave that propagates backward. The non-local constructive interaction of such waves gives rise to a self-sustained global in time instability, which, in some circumstances, is able to radiate an intense acoustic field.

St (n) β = ω (n) β 2π (D 1) (a) (b)
In the following analysis, we reconsider the discharge coefficient coefficient α vena = ρU 2 M 2(Pin-Pout) , which can be thought as a measure of the vena contracta phenomenon, assuming that the jet contracts to a top-hat jet with constant velocity U J and radius R J . Then applying Bernouilli law, one obtains

α vena = U M U J = πR 2 J πR 2
h that was introduced in [START_REF] Fabre | The acoustic impedance of a laminar viscous jet through a thin circular aperture[END_REF] in the discussion of the work of [START_REF] Howe | On the theory of unsteady high reynolds number flow through a circular aperture[END_REF], which is classically associated with the pressure loss across the aperture, and it relates the mean velocity U M with the jet velocity U J . In Fabre et al. (2019, Sec. 2.5 & App. A) it is argued the importance of the vena contracta coefficient α vena and the actual value of the phase velocity of the Kelvin Helmholtz instability, which depends on the frequency ω but which is around U c = U J /2 for sufficiently high ω. The value of the discharge coefficient tends asymptotically to α vena → 0.61 for large Re, however in the range of Reynolds numbers where transition occurs α vena is maximal and takes values α vena ∈ [0.7, 0.76]. In the following, we consider a constant discharge coefficient α vena = 0.76, which was the value reported in Fabre et al. (2020, Fig. 6) for most of β in the range of Reynolds numbers where the transition occurs.

In order to estimate the frequency shift ∆St β let consider the travel time of each travelling wave past the hole. The hydrodynamic travelling wave takes ∆τ KH = L h Uc =

2αvenaL h U M
and the acoustic wave ∆τ ac = M L h U M : therefore the total travelling time

∆τ = ∆τ KH + ∆τ ac = L h (2α vena + M ) U M ≈ 2α vena L h U M , f (n) = 2n -1 2∆τ , ∆f = 1 ∆τ (D 2)
where it is considered that M ≪ 1 and the convective velocity of the hydrodynamic perturbation is 2αvena U M as it is displayed in fig. 19. Thus, the associated Strouhal shift

∆St β = ∆f L h U M ≈ 1 2α vena . (D 3)
Similarly, the Strouhal frequency of each C n branch is estimated as

St (n) β = 2n -1 4α vena , n = 1, 2, 3 . . . ( D 4)
In the previous reasoning, it has been implicitly assumed that only odd mode structures as those depicted in Fabre et al. (2020, Fig. 9) lead to a conditional instability. The superposition of the base flow with odd modes (resp. even) at the instant of the cycle where the flow rate through the hole is maximum results in an upward (resp. downward) displacement of the shear layer, thus pressure is increased (resp. decreased) in the presence of odd modes (resp. even). This implies that the hole impedance is negative, which is the criterion of existence of conditional stability. Finally, let us discuss the influence of compressibility regarding modifications in frequency at the threshold of instability for each of the unstable modes. Compressibility induces a weak variation in the vena contracta coefficient, see fig. 16 (a), while the threshold of instability is significantly modified for a given value of β near the left asymptote for conditional instabilities (fig. 18 (a)) and hydrodynamic instabilities (fig. 11 (a)). Thus, for short holes β (here short refers to a value of β near a vertical asymptote of the corresponding instability), the critical Reynolds number is considerably modified by compressible effects, which in turn induces a change in the vena contracta coefficient (more important than the variation of the vena contracta by compressibility at constant Reynolds number). For those cases, one could evaluate the shift in frequency δ| where the first term, the variation in the vena contracta due to a variation in Mach, is negligible with respect to the variation in Mach number for long holes, and thus, the variation of frequency for those instabilities is of lesser importance.

Zmax Rmax Z0 Lc γz,c R0 Rc γr,c ntri M1 50 50 10 2.5 0.3 5 2.5 0.3 42460 tri. M2 300 300 10 5 0.3 10 2.5 0.3 165841 tri. M3 600 600 10 5 0.3 10 2.5 0.3 235874 tri. M4 600 600 170 5 0.3 170 5 0.3 314525 tri. M5 1500 1500 10 5 0.3 10 2.5 0.3 413356 tri.

Table 1: Meshes used for cases where the complex mapping technique is adopted. Note: the smooth transition functions are defined as g z (Z) = tanh Z-Z0 Lc

Appendix E. Computational domains & and absorbing boundary layer

This section discusses the design of computational domains used for the computation of steady-states and eigenmodes with the full compressible formulation, and steadystates and the forced harmonic response with the incompressible formulation. The computational strategy must be designed in such a way to avoid the presence of spurious eigenvalues/eigenmodes. The creation of meshes for the full compressible formulation follows a block-structured strategy, similar to the one sketched in fig. 4. The domain is divided into three regions: an inner region with the highest vertex density, a mid region with intermediate vertex density, and a coarser region for the absorbing layer. Meshes employed for the incompressible case are composed of two regions, a physical domain with the highest vertex density and a coarser region for the absorbing layer. Table 1 lists a number of computational meshes employed in this study (M 2 , M 3 , M 4 and M 5 for the full compressible case, and M 1 for the incompressible case). However, such a list is not exhaustive because in addition to the block-strategy refinement, the computational domain has been locally refined following an adaptive local refinement (AMR) procedure, where the metric for the refinement is based on the steady-state and on the eigenmodes (resp. steady-state and forced harmonic response in the incompressible case), cf Hecht (2012); [START_REF] Fabre | A practical review on linear and nonlinear global approaches to flow instabilities[END_REF].

The absorbing boundary layer corresponds to the complex mapping technique [START_REF] Sierra | Efficient stability analysis of fluid flows using complex mapping techniques[END_REF] where a coordinate transformation G is defined as follows:

G z : R → C such that z = G z (Z) = 1 + iγ z,c g z (Z) Z G r : R → C such that r = G r (R) = 1 + iγ r,c g r (R) R (E 1)
Here g z (Z) (resp. g r (R)) has to be chosen as a smooth function such as g z (Z) = 0 for Z < Z 0 and g z (Z) ≈ 1 for Z > Z 0 + L c up to Z max . The complex mapping acts on a finite region of length z CM = Z max -Z 0 + L c . The function g z is defined as

g z (Z) = tanh Z-Z0 Lc 2
. The application of this map to the linearized Navier-Stokes requires that each spatial derivative, within the complex mapped region, is modified as follows:

∂ ∂z ≡ H z ∂ ∂Z with H z (Z) = ∂G z ∂Z -1
.

(E 2)

An example of parameters for the usage of the complex mapping layer for the incompressible case is listed in table 1 (M 1 ). In the full compressible formulation, one must pay particular attention about the extension of the complex mapping region, which here is selected to cover at least two acoustic wavelengths, i.e. z CM = r CM > 1 StM . For instance, the largest acoustic wavelength in this study corresponds to the validation case in section 5.1, where λ ac ≈ 600 (St ≈ 1/3 and M = 5 • 10 -3 ), for which M 5 is an appropiate choice.

Figure 1 :

 1 Figure 1: Sketch of the open/open configuration.

Figure 2 :

 2 Figure 2: Sketch of the cavity/open configuration.

  (a): For the open/open configuration, Z h = -2Z rad , or equivalently :

Figure 3 :

 3 Figure 3: Linear scale representation of the zero computation of eq. (3.5) with the (a) zeroth-order eq. (3.7) and the (b) first order eq. (3.8) approximations.

Figure 4 :

 4 Figure 4: Schematic representation of the computational mesh for both configurations, (a) open/open case, (b) closed/open case: z -∞ , z ∞ , r ∞ are, respectively, the location of the physical inlet, outlet and lateral boundaries. The physical domain is padded into a complex mapping layer with a radial extension r CM (resp. axial z CM extension). The inner domain corresponds to an inner region with the highest vertex density: z -, z + , r + are, respectively, the location of the left, right and lateral boundaries of this inner domain; in the closed/open case the inner domain includes the cavity located upstream of the hole.

  (a), composed of three regions: an inner domain with the highest vertex density, the physical domain and an absorbing layer to eliminate the appearance of spurious eigenvalues. The absorbing layer corresponds to the complex mapping technique, cf. Sierra et al. (2020). The boundary conditions of the linearized full compressible formulation for the open/open case are as follows

Figure 5 :

 5 Figure 5: Growth rate (left) and frequency (right) of eigenmodes as a function of V in , M and Re for β = 0.3. Lines were obtained from the matched asymptotic model and points with the compressible LNSE. Solid lines denote unstable regions, dashed lines are used for stable zones.

Figure 6 :

 6 Figure 6: C 1 eigenmode for Re = 1200 (real part in upper region and imaginary part in the lower region) at M = 5 • 10 -3 . (a) Pressure, (b) temperature and (c) density.

Figure 7 :

 7 Figure 7: Real part of the pressure component of higher order cavity modes for M = 2 • 10 -1 and a cavity of V in = 10 4 .

Figure 8 :

 8 Figure 8: Regions of conditional stability in the (V in , M ) plane for β = 0.3.

Figure 9 :

 9 Figure 9: Regions of conditional stability in the (V in , M ) plane for β = 1.

Figure 10 :

 10 Figure 10: Regions of conditional instability in the (χ, Re) plane.

Figure 11 :

 11 Figure 11: (a) Neutral curves of instability of hydrodynamic modes H 2 (solid lines) and H 3 (dashed lines) for M = 10 -2 , 5 • 10 -2 and M = 10 -1 . (b) Strouhal evolution of the H 2 modes with β at the threshold of the instability. (c) The same as (b) for the H 3 modes. Dots correspond to incompressible results of Fabre et al. (2020), also reported in appendix D.

Figure 12 :

 12 Figure 12: Sensitivity to Mach number variations of the complex frequency ω at the threshold of instability. Dashed line denotes the impedance estimation. (a) Mode H 2 for β = 2. (b) Mode H 3 for β = 0.8.

Figure 13 :

 13 Figure 13: Real part of the pressure component of neutral eigenvalues H 2 (upper part) and H 3 (lower part) at M = 5 • 10 -2 . (a) β = 1 and (b) β = 2.

  Figure13(a) displays the pressure levels in logarithmic scale for β = 1. In that figure one can appreciate a monopolar-like radiation for H 3 , however at M = 5 • 10 -2 the H 2 mode displays a radiation pattern with a preferential direction aligned with the stream-wise coordinate. In such a configuration, either the hole and the jet emit sound, which produces an uneven radiation of sound downstream and upstream of the hole. These observations can be better appreciated in fig.14 (a), where it is clear the departure from a isotropic radiation for M = 10 -2 and the sound emission for M = 5 • 10 -2 . For β = 2 ( fig.13(b) and fig.14 (c) ) the neutral eigenmodes display a fairly monopolar-like radiation; for the H 3 mode (fig.14 (d)), higher pressure levels are measured downstream of the hole for directions forming an angle less than 45 degrees with the axis of symmetry.

Figure 14 :

 14 Figure 14: Directivity patterns, in logarithmic scale, measured at r x = √ r 2 + x 2 = 150. (a) H 2 and β = 1; (b) H 3 and β = 1; (c) H 2 and β = 2; (d) H 3 and β = 2. Colour legend: M = 10 -2 , M = 2 • 10 -2 and M = 5 • 10 -2 .

Figure 15 :

 15 Figure 15: Contour plot of the baseflow q 0 for β = 0.3, Re = 1400 and M = 5 • 10 -2 . (a) Spatial evolution of the Mach number. (b) Spatial evolution of the sensitivity of the axial velocity with respect to the Mach number ∇ M w 0 .

Figure 16 :

 16 Figure 16: Discharge coefficient as a function of Mach number for several Re for β = 0.3. Dots correspond to computed values of α vena , lines correspond to the theoretical estimation, cf. (Bragg 1960). (b) Axial velocity profile of the sensitivity with respect to Mach number for β = 0.3, Re = 1600, M = 10 -1 at z = 0.2.

Figure 17 :

 17 Figure 17: (a) The complex contour Γ of integration enclosing the unstable complex plane, where Cauchy's argument principle is applied. (b) Nyquist curve for the augmented system impedance along Γ .

  the first-order Taylor expansion at ω 0 provides 0

Figure 18 :

 18 Figure 18: a) Thresholds of the conditional instabilities C 1 , C 2 , C 3 and C 4 as a function of β and Re. b) Strouhal evolution with β. Solid lines are the boundaries of the conditional stability computed with incompressible LNSE and dashed lines are the estimation with vena = 0.76

  M ∆St β = ∆St β (M + ∆ M ) -∆St β (M ) as δ| M ∆St β = -2 α vena (Re c (M + ∆M ), M + ∆M ) -α vena (Re c (M ), M ) + ∆M 2α vena (Re c (M ), M ) + M 2 (D 5)

.

  In the followingL max = z -∞ + z CM = z ∞ + z CM and R max = r ∞ + r CM (fig.4).

† Email address for correspondence: fgiannetti@unisa.it

Appendix A. Computation of the hole impedance thanks to

Incompressible linearized Navier Stokes equation

A.1. Incompressible Navier-Stokes equations

Under the hypothesis of acoustic compactness, discussed in section 2.3, the flow is assumed to be locally incompressible in the region of the hole, where the fluid motion is governed by the incompressible Navier-Stokes equations:

The stability of the steady-state u 0 , p 0 is investigated by using the linearized approach, in which the total flow field is decomposed into the sum of a steady base flow and a small time-harmonic perturbation as:

. Incompressible Navier-Stokes equations -Base-flow equations

The base flow is the solution of the steady version of the Navier-Stokes equations:

N S[u 0 ; p 0 ] = 0 with the following set of boundary conditions:

This problem is solved using a classical Newton iteration.

A.3. linearized incompressible Navier-Stokes equations -Forced problem

The linear perturbation is governed by the following equations:

where LN S 0 is the linearized Navier-Stokes operator around the base flow and B is a weight operator defined as follows:

Equation (A 5) is complemented with the following boundary conditions :

A nonzero perturbation of the flow rate q ′ is imposed, fixed arbitrarily, to q ′ = 1. A 6 thus leads to a non-homogeneous Dirichlet boundary condition at the inlet plane, treated by imposing a constant axial velocity ûx . The problem can be symbolically written as [LN S + iωB] [û; p] = F, (A 8)