- Email: [email protected]
The propagation of acoustic waves in a slowly varying duct with radially sheared axial mean flow
Journal of Sound and Vibration 332 (2013) 3937–3946
Contents lists available at SciVerse ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
The propagation of acoustic waves in a slowly varying duct with radially sheared axial mean flow A.E.D. Lloyd, N. Peake n Department of Applied Mathematics & Theoretical Physics, University of Cambridge, Wilberforce Road, Cambridge CB3 0WA, United Kingdom
a r t i c l e i n f o
abstract
Article history: Accepted 1 February 2013 Handling Editor: P. Joseph
We consider the propagation of acoustic waves along a cylindrical duct carrying radially sheared axial mean flow, in which the duct radius is allowed to vary slowly along the axis. In previous work [A.J. Cooper & N. Peake, Journal of Fluid Mechanics 445 (2001) 207–234.] radially sheared axial mean flow with nonzero swirl in a slowly varying duct was considered, but in this paper we set the swirl to zero, thereby allowing simplification of the calculations of both the mean and unsteady flows. In this approach the acoustic wavenumber and corresponding eigenfunction are determined locally, while the wave amplitude is found by solving an evolution equation along the duct. Sample results are presented, including a case in which, perhaps surprisingly, the number of cut-on modes increases as the duct radius decreases. & 2013 Elsevier Ltd. All rights reserved.
1. Introduction There has been considerable interest in the propagation of acoustic waves in flow-bearing ducts for many years, motivated by a number of different industrial applications. In the context of turbomachinery noise, one is often interested in the two different situations; one situation in which the mean flow is both radially sheared and contains a significant swirling component, as found between a rotor row and a stator row; and a second situation in which the swirl is negligible but the axial flow is still radially sheared, as found in the bypass duct downstream of the outlet guide vanes. In both the cases the presence of nonzero mean vorticity is significant and must be included in the modelling. Also in both the cases, but especially in the latter case, the axial variation of the duct cross-section should also be considered. The propagation of acoustic waves in a straight duct carrying axial flow has received much attention, including the seminal work of Pridmore-Brown [1]. Many other early references are given in [2]. A common issue throughout has been the role of wall linings in attenuating the sound, and we would like to mention here in particular an ingenious early paper by Doak and Vaidya [3], in which an approximate analytical expression for attenuation coefficients, valid in a range of geometries, was derived. More recently, Vilenski and Rienstra [4,5] have completed a detailed numerical study over a wide range of frequencies. Although acoustic (pressure-dominated) waves are of primary importance, the sheared axial flow also allows the propagation of hydrodynamic (vorticity-dominated) waves, described in inviscid flow by the critical layer. The inclusion of mean swirl has been considered, for instance by [6–9]; swirl has a significant quantitative effect on the acoustic field, and also gives rise to an additional family of hydrodynamic, ‘nearly convected’ modes located close to the critical layer in wavenumber space. Full details of these hydrodynamic features are given in the above references, and a recent review is provided in [10]. Axially varying ducts have been considered, and in particular the work of Rienstra [11] on slowly varying ducts with radially uniform axial mean flow has proved very influential - see [12–17] for a series of extensions. n
Corresponding author. Tel.: þ 44 1223 339058. E-mail address: [email protected] (N. Peake).
0022-460X/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.jsv.2013.02.038
3938
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
r
r=R2 r=R1
θ X
Fig. 1. Schematic of slowly varying duct.
The particular extension of [11] that we wish to highlight here, however, is the case of acoustic wave propagation in a slowly varying duct containing both mean swirl and radially sheared axial mean flow, which was completed by [18] (with ramifications for noise generation presented in [19]). The analysis in [18] is necessarily complicated, due to the rather general nature of the vortical mean flow considered. The bypass-duct situation of zero swirl but radially sheared axial flow is certainly contained as a special case of the theory in [18], but is not explored there. In fact, by setting the swirl to zero in [18] it turns out that the determination of both the mean and the unsteady flow is simplified, and this is the subject of the current paper. In Section 2 we show how the evolving mean flow can be calculated along the duct, and present the unsteady flow equations; in Section 3 we show how to determine the evolution of the acoustic modes along the slowly varying duct, and in Section 4 we give sample results. In many situations it can be expected that the number of cut-on modes decreases as the duct radius decreases (see for instance [20]), but we present a case in which the number of cut-on modes increases as the duct radius decreases. 2. Problem formulation Consider an annular duct with circular cross section, described using a cylindrical coordinate system ðx,r,θÞ, as shown in Fig. 1. The duct carries a compressible, perfect, isentropic flow, which is made up of a steady mean component plus a small unsteady perturbation. Throughout this paper we non-dimensionalise all lengths by the outer duct radius at x¼0, and velocities and densities by the mean sound speed and mean density on the outer wall at x ¼0, respectively. The duct crosssection varies slowly in the axial direction and we introduce a slow axial scale, X ¼ ϵx, where ϵ⪡1 is a small parameter defined as a measure of the axial slope of the duct walls. The duct inner and outer radii are denoted by R1,2 ðXÞ. The flow has total velocity v, pressure p, density ρ and sound speed c, and is decomposed in the form ~ c~ , ρ, ~ ~ pÞðX,r,θ,tÞ, ðv,c,ρ,pÞðX,r,θ,tÞ ¼ ðV,C,D,PÞðX,rÞ þðv,
(1)
~ c~ , ρ, ~ is the small-amplitude (linear) unsteady perturbation. We will describe how ~ pÞ where ðV,C,D,PÞ is the steady flow and ðv, to determine these two components in the next two subsections. 2.1. Mean flow We consider a mean flow field which has zero swirl but which has axial and radial components which vary with both X and r. It then follows that the mean flow has the asymptotic form ðU,V,C,D,PÞðX,r; ϵÞ ¼ ðU 0 ðX,rÞ,ϵV 1 ðX,rÞ,C 0 ðX,rÞ,D0 ðX,rÞ,P 0 ðX,rÞÞ þ Oðϵ2 Þ,
(2)
while the assumption of a homentropic perfect gas relates P0 and C0 to D0 via P0 ¼
1 γ D γ 0
ðγ−1Þ=2
C 0 ¼ D0
and
:
(3)
At OðϵÞ the continuity equation is ∂ðD0 U 0 Þ 1 ∂ðrD0 V 1 Þ þ ¼ 0, ∂X r ∂r
(4)
and can be satisfied by writing velocity components U0 and V1 in terms of the streamfunction ψðX,rÞ, such that U0 ¼
1 ∂ψ , rD0 ∂r
V1 ¼ −
1 ∂ψ : rD0 ∂X
(5)
Bernoulli's relation (enthalpy is constant on streamlines) can now be written as γ−1
1 2 D0 U þ ¼ HðψÞ, 2 0 γ−1
(6)
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3939
where the enthalpy H is a function of ψ specified by the conditions upstream. For homentropic flow Crocco's relation is
where the mean vorticity ξ is
and taking the r component of (7) gives
V ξ ¼ ∇H,
(7)
∂U 0 ∂V 1 þ ϵ2 eθ , ξ¼ − ∂r ∂X
(8)
∂ 1 ∂ψ ¼ rD0 H′ðψÞ: ∂r rD0 ∂r
(9)
We note here that our non-swirling mean flow Eqs. (6) and (9) take a simpler form than their nonzero-swirl counterparts, Eqs. (17) and (21) in [18]. Now differentiating (6) with respect to r and substituting for H′ðψÞ from (9) we see that ∂D0 =∂r ¼ 0, and therefore that the mean density, and hence from (3) also the mean pressure and sound speed, are functions of X only. The boundaries of the duct are streamlines of the mean flow, along which ψðX,rÞ is constant; the values of these constants will be specified as a part of the initial conditions at X ¼0. For convenience the axial velocity profile at X¼ 0 is written as U 0 ðX ¼ 0,rÞ ¼ g′ðrÞ=r, and the function g(r) is chosen to give a suitable initial profile. Using this profile and integrating Eq. (5) we find ψðX ¼ 0,rÞ ¼ gðrÞ, with the constant of integration taken to be zero without loss of generality (note that our nondimensionalisation implies that D0 ð0Þ ¼ 1). We make the assumption that there is a streamline connecting each point downstream with a point at X¼0 (i.e. no regions of reverse mean flow), allowing us to equate the Bernoulli quantity given by Eq. (6) at X¼0 to the value on the same streamline at an arbitrary distance downstream. After rearrangement this gives rise to the equation ! 1=2 ∂ψ g′ðr s Þ2 2 γ−1 ¼ D0 ðXÞr 1−D0 ðXÞ þ , (10) ∂r γ−1 r 2s where the radial coordinate r s ðX,ψÞ is such that the streamline passing through the point (r,X) also passes through the point ðr s ,0Þ. We now integrate (10) to give " # Z ψ g′ðr Þ2 −1=2 2D20 r 2 −R21 ðXÞ s : (11) þ 1−Dγ−1 dψ ¼ 0 r 2 γ−1 s gðR1 ð0ÞÞ We first set r ¼ R2 ðXÞ in (11) (so that ψ-gðR2 ð0ÞÞ), giving an implicit equation for D0 ðXÞ which can be solved iteratively by marching downstream from X¼ 0. Once D0 ðXÞ has been determined, ψðr,XÞ can be found iteratively from (11), and then finally the mean flow determined from (5). Note that in [18] the calculation of the mean flow with nonzero swirl is more complicated than this: in nonzero swirl the mean pressure, density and sound speed depend on r as well as X; determination of the streamfunction as a function of position requires solution of two coupled ordinary differential equations; and while in principle the method presented in [18] is valid for general upstream conditions, for simplicity only a straightforward case of solid-body rotation and uniform mean flow upstream were considered there. As specific examples we consider two annular ducts, one slowly contracting and the other slowly expanding, based on the examples used in [18]. For the contracting case the duct boundaries are defined by R1 ðXÞ ¼ 0:4482 þ0:05 tanhð2X−2Þ,
R2 ðXÞ ¼ 0:9518−0:05 tanhð2X−2Þ,
(12)
and for the expanding case the duct boundaries are defined by R1 ðXÞ ¼ 0:45180−0:05 tanhð2X−2Þ,
R2 ðXÞ ¼ 1:04820 þ0:05 tanhð2X−2Þ:
(13)
We note that the duct mean line, r ¼ ðR1 ðXÞ þ R2 ðXÞÞ=2, is straight in both the cases. Clearly, other examples are possible in which duct area variation is combined with curvature of the duct mean line, but these additional cases will not be considered here. Following [4] we model the axial velocity profile at X¼ 0 as a parabola, chosen so that the maximum velocity lies in the centre of the duct half way between the inner and outer walls U 0 ðX ¼ 0,rÞ ¼ U max
1 þ sðr−R1 ð0ÞÞð1−rÞ 1 þsð1−R1 ð0ÞÞ2 =4
,
(14)
and where U max ¼ 0:5 and s ¼ 8, corresponding to U 0 ðR1 ð0ÞÞ ¼ 0:76U max . This choice of U0 corresponds to gðrÞ ¼ U max
6r 2 ð1−R1 ð0ÞsÞ þ4sð1 þ R1 ð0ÞÞr 3 −3sr 4 12 þ 3sð1−R1 ð0ÞÞ2
:
(15)
In Fig. 2 we show the mean flow for the contracting duct described in Eq. (12). Note that V1 is small, and therefore that the radial velocity, ϵV 1 , is of negligible magnitude in comparison to the axial component U0. However, the figure also shows that U0 varies significantly along the contracting duct, in which a significant acceleration downstream is accompanied by a reduction in D0.
3940
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
Fig. 2. Leading order components of (a) the axial mean velocity, U 0 ðX,rÞ; (b) the radial mean velocity, V 1 ðX,rÞ. The bold lines denote the duct walls. The axial variation of (c) the mean density, D0 ðXÞ; (d) the mean sound speed C 0 ðXÞ. All for a contracting duct with an initial parabolic mean flow profile (from the left). In (a) the contours are plotted at intervals of 0.05 and in (b) at intervals of 0.005.
2.2. Unsteady flow equations Having obtained the solution for the mean flow we now consider the unsteady disturbance field. Following the work of [7] we decompose the perturbation velocity into a sum of vortical and potential components v~ ¼ uR þ∇ϕ:
(16)
This decomposition is made unique by requiring that the unsteady pressure is expressed solely in terms of the unsteady potential p~ ¼ −D0
Dϕ , Dt
(17)
where D=Dt is the convective derivative D=Dt ¼ ∂=∂t þV⋅∇. The Euler and continuity equations now become DuR þðuR ∇ÞV ¼ −ξ ∇ϕ, Dt ! D 1 D 1 1 − ∇ ðD0 ∇Þ ϕ ¼ ∇ D 0 uR , Dt C 20 Dt D0 D0
(18)
(19)
where ξ is the mean vorticity as defined in (8). In the absence of mean vorticity equation (18) decouples from (19), and uR can be determined everywhere by integrating (19) downstream starting from given conditions at upstream infinity. However, in the present case the unsteady vorticity is fully coupled to the unsteady pressure through the nonzero mean vorticity, and Eqs. (18) and (19) must be solved simultaneously. We will impose boundary conditions to model the presence of an acoustic liner on the inner and outer duct walls of given impedance Z1 and Z2 respectively. For an arbitrary curved wall the correct boundary condition was first given by [21], and has been implemented by [11,18]. It is p~ iω v~ nj ¼ iωþ V ∇−nj ⋅ðnj ⋅∇VÞ , (20) Zj
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3941
where nj , j ¼ 1,2 are the outward (i.e. pointing out of the fluid) unit normal vectors at the wall. For the slowly varying duct these normals are n1 ¼ −
n2 ¼
er −ϵR′1 ex ð1 þϵ2 R′21 Þ1=2 er −ϵR′2 ex
ð1 þ ϵ2 R′22 Þ1=2
,
:
(21)
(22)
Hard-wall boundary conditions are given by the limit Z j -∞. The unsteady flow equations permit two essentially distinct families of solutions, corresponding to acoustic modes which are pressure-dominated and hydrodynamic motions which travel close to the local mean flow speed. In this paper we consider only the acoustic modes. 3. Acoustic modes As in [18] we use the method of multiple scales and seek a solution of the form ðϕ,uRx ,uRr ,uRθ Þðx,r,θ,t; ϵÞ ¼ ðA,H,R,TÞðX,r; ϵÞE:
(23)
Here Z i E ¼ exp iωt−imθ− ϵ
X
kðηÞ dη ,
(24)
ω is the frequency of the acoustic perturbation, m is the azimuthal wavenumber and k(X) is the local axial wavenumber. The amplitudes and the axial wavenumber evolve along the slowly varying duct and therefore depend on the slowly varying X coordinate. Substituting Eq. (23) into Eqs. (18) and (19) gives rise to four coupled equations, which, up to and including OðϵÞ terms, become ! " ! !# ! ∂2 A 1 ∂A Λ2 m2 2 ∂R imT ϵ ∂ iϵ ∂ U0Λ 1 ∂ rD0 V 1 ΛA2 2 , þR þ − þ þ − −k Aþ ðD0 H Þ−ikH ¼ þk D0 A þ r ∂r ∂r r D0 ∂X D0 A ∂X r ∂r ∂r 2 C 20 r 2 C 20 C 20 (25) ∂U 0 ∂R ∂V 1 ∂R ∂A ∂U 0 þR − ¼ −ϵ V 1 þU 0 , iΛRþ ikA ∂r ∂X ∂X ∂r ∂r ∂r
(26)
∂T V 1 T ∂T þ þ U0 , iΛT ¼ −ϵ V 1 ∂r r ∂X
(27)
∂A ∂U 0 ∂U 0 ∂H ∂H þU 0 , ¼ −ϵ H þV1 iΛH þ R þ ∂r ∂r ∂r ∂X ∂X
(28)
where Λ ¼ ω−kU 0 . Similarly, the boundary conditions (20) give ! dRj ∂A D0 Λ2 A iϵ ∂ ∂ ∂V 1 dRj ∂U 0 D0 ΛA2 þR 7 U0 þV1 − ðkA þiHÞ 7 þ iω ¼ ϵω ∂r Zj A ∂X ∂r ∂r dX dX ∂r Zj
(29)
on r ¼ Rj ðXÞ, j¼1,2, where 7 refers to R1 and R2. Eqs. (25)–(29) agree with the corresponding Eqs. (36)–(40) in [18] when the mean swirl is set to zero. We now proceed to expand the amplitudes in (23) in powers of ϵ ðA,H,R,TÞðX,r; ϵÞ ¼ ðA0 ,H 0 ,R0 ,T 0 ÞðX,rÞ þϵðA1 ,H 1 ,R1 ,T 1 ÞðX,rÞ þ Oðϵ2 Þ,
(30)
and substituting into Eqs. (36)–(39) we find at Oð1Þ and OðϵÞ Lψ 0 ¼ 0,
(31)
Lψ 1 ¼ f,
(32)
respectively, where L ¼ L1 −kL2 , and the eigenfunctions take the form ψ n ¼ ðAn ,ηn ,Rn ,iH n Þ for n¼0,1, with ηn ¼ kβ20 An and β20 ¼ 1−U 20 =C 20 . Explicit expressions for L1 , L2 and f are given in Appendix A. We note that the azimuthal component of the perturbation vortical velocity is always equal to zero (i.e. T 0 ¼ T 1 ¼ 0); in the case of nonzero mean swirl [18] all the three components of the vortical velocity are nonzero and must be determined. Eqs. (31) and (32) are to be solved subject to the boundary conditions ∂A0 D0 ðω−kU 0 Þ2 iω þR0 7 A0 ¼ 0 (33) ∂r Zj
3942
and
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
! dRj ∂A1 D0 ðω−kU 0 Þ2 A1 i ∂ ∂ ∂V 1 dRj ∂U 0 D0 ðω−kU 0 ÞA20 þV1 − þ R1 7 ðkA0 þiH 0 Þ 7 þ ¼ω U0 iω A0 ∂X ∂r ∂r ∂r Zj dX dX ∂r Zj
(34)
respectively, which are obtained at OðϵÞ and Oðϵ2 Þ from (29). The leading-order problem defined by (31) and (33) is solved for the local axial wavenumber k(X) as an eigenvalue, with corresponding eigenfunction ψ^ 0 ðX,rÞ. At this order the solution is therefore known up to a slowly varying amplitude function, N(X), so that ψ 0 ðX,rÞ ¼ NðXÞψ^ 0 ðX,rÞ:
(35)
To determine N(X) we will take the inner product of (32) with the homogeneous solution of the adjoint problem. The inner product between the vectors A ¼ ðA1 ,A2 ,A3 ,A4 Þ and B ¼ ðB1 ,B2 ,B3 ,B4 Þ is defined as Z R2 4 ∑ Ann Bn r dr, (36) 〈A,B〉 ¼ R1
n¼1
and we denote the solution of the adjoint problem as ψ †0 ¼ ðY 1 ,Y 2 ,Y 3 ,Y 4 Þ. With this choice of inner product the adjoint problem is found to be L† ψ †0 ¼ 0,
(37)
†
where the adjoint operator L is defined in Appendix A, subject to the boundary condition " # ∂Y n1 D0 ðω−kU 0 Þ2 Y n1 n ∂U 0 Y4 ¼ 0: −D0 7 ∂r ∂r iωZ j
(38)
We therefore have the two systems of equations to solve, one for the eigenvector given in (31) and (33), the other for its adjoint given in (37) and (38). By comparing the two problems, we find that the adjoint eigenvector ψ †0 can be related to the direct eigenvector ψ 0 by ! 2ωU 0 iD0 Hn D0 Rn0 n : (39) þ k An0 , Y 3 ¼ − 0 , Y 4 ¼ − Y 1 ¼ An0 , Y 2 ¼ D0 β20 2 2 ∂U ∂U 0 0 C 0 β0 ∂r ∂r The inner product of (32) with ψ †0 gives 〈ψ †0 ,f〉 ¼ 0, and after some algebra this leads to the governing equation for N(X) in the form FðXÞ
d ðN2 ðXÞÞ ¼ N2 ðXÞGðXÞ: dX
Here F(X) and G(X) are defined as Z FðXÞ ¼ and
R2 ðXÞ
r R1 ðXÞ
U0 Λ C 20
2
! þ k D0 A^ 0
(40)
2
3 2 3 2 2 ^ 02U0 ^ 02U0 D R Λ A D R Λ A 2 1 5 þ4 0 5 dr þ 4 0 iωZ 2 iωZ 1 R2
# n n ∂A^ 0 ^ ^ ^ ^ ^ þ iY 3 f 3 r þ iY 4 f 4 r dr r irD0 H 0 ∂X C 20 R1 ðXÞ R1 ðXÞ 2 0 13 2 0 13 ^2 ^2 D R ∂ ∂ ∂V dR ∂U D Λ A D R ∂ ∂ ∂V dR ∂U D Λ A 0 2 1 2 0 0 0 1 1 1 0 0 0 0 @ A 5 −4 @ A5 þV1 − þV1 − U0 þ U0 þ −4 ∂X ∂r ∂r ∂X ∂r ∂r iω dX ∂r Z2 iω dX ∂r Z1 d GðXÞ ¼ − dX
Z
R2 ðXÞ
"
U0 Λ
!
(41) R1
# Z 2 ^ ^ ^ þ k D0 A 0 þ irD0 A 0 H 0 dr þ
R2
R2 ðXÞ
"
(42) R2
respectively. The functions f^ 3 and f^ 4 refer to the components of the vector f given in Appendix A, but with non-hatted variables replaced with their hatted equivalents. The general solution to (40) is now simply Z X GðηÞ dη , (43) N2 ðXÞ ¼ N 2 ð0Þ exp 0 FðηÞ where we have chosen to fix the amplitude relative to its value at X¼0. Eq. (43) provides the most convenient form for evaluation of the variation of N(X). However, for ease of comparison with the radially uniform mean flow solution in [11], we also note that we can rearrange (40) into the form #
Z R2 ðXÞ Z R2 ðXÞ " n n d ∂A^ 0 2 2 d 2 ^ ^ ^ ^ ^ ^ ^ ðFN Þ ¼ −N þ iY 3 f 3 r þ iY 4 f 4 r dr: iD0 A 0 H 0 dr þN irD0 H 0 (44) dX dX ∂X R1 ðXÞ R1 ðXÞ In the special case of potential disturbances in radially uniform mean flow, U 0 ðX,rÞ is replaced by U 0 ðXÞ and the vortical component of the unsteady velocity becomes zero (when the mean vorticity is zero, zero unsteady vorticity upstream
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3943
implies zero unsteady vorticity throughout the flow). The right hand side of Eq. (44) then becomes exactly zero, and (44) agrees exactly with the corresponding equation (4.10) in [11]. 4. Results Our first step is to solve numerically the leading-order eigenvalue problem L1 ¼ kL2 , together with the boundary condition (33), to obtain the unknown axial wavenumbers k and the corresponding eigenfunctions, at each X. As in [18], to do this we employ a Chebyshev spectral collocation method of the form given by [22], in which we choose nr collocation points in the radial direction. Full details of this method can be found in these references and in [23]. One of the well-known characteristics of pseudo-spectral methods is the appearance of spurious eigenvalues, for example see [24]. In our problem these include spurious modes located in an egg-shaped region around the critical layer, as noticed by [25]. Following the work of [14], we apply two filtering procedures to remove these spurious eigenvalues; the first removes eigenvalues whose position varies significantly as nr is varied, while the second removes eigenvalues for which the corresponding eigenfunction is not sufficiently well resolved. Full details are again given in [23]. We also note that if one were interested in identifying hydrodynamic modes, as opposed to the acoustic modes considered throughout this paper, then even more care would have to be taken when separating genuine and spurious eigenvalues. For instance, in work on swirling ducted flow [9], classical numerical shooting was found to be effective, even for modes close to the critical layer. This technique could be included as a further check on the results if the above pseudospectral methods were applied to the hydrodynamic modes. We now look in detail at the axial variation of the spectra for the expanding duct (Fig. 3) and the contracting duct (Fig. 4). In both the cases, there are two cut-on modes propagating downstream and two cut-on modes propagating upstream at X¼0. As we move downstream the cut-on and cut-off eigenvalues move in the complex k plane. In the expanding case there are still just two upstream and two downstream propagating cut-on modes at X ¼1.2, and the modes which were cut-off at X¼0 remain cut-off at X¼ 1.2—in fact this continues to be the case as X-∞. However, more complicated behaviour is apparent in the contracting case, in which the less cut-off modes at X¼0 move significantly towards the real axis downstream, and indeed by X¼1.6 (Fig. 4) two pairs of original cut-off modes have moved onto the real axis, so that there are now four cut-on modes propagating downstream and four cut-on modes propagating upstream. This is perhaps surprising, since one normally thinks of modes becoming cut-off as they propagate along a narrowing duct; see for instance [20], in which a model for acoustic resonance in an intake duct involves trapping between the duct narrowing upstream and the swirl behind the fan downstream. To understand this behaviour, we consider the special case in which the duct is hollow (i.e. R1 ðXÞ ¼ 0) and the flow far upstream is uniform (i.e. s ¼ 0 and U 0 ðX,rÞ ¼ U 0 ðXÞ, as studied by [11]). In this case there are well-known analytic expressions for the axial eigenvalues (see for example [11, Eq. 4.9]), and if we simplify further to the case of hard walls we find that the number of cut-on modes corresponds to the number of zeros of J′m ðzÞ, j′mn for n ¼ 1,2,…, which satisfy j′mn o
ωR2 C 0 ð1−U 20 =C 20 Þ
:
(45)
In (45) the X variation of C 0 ðXÞ is rather weak (see Fig. 2b) and can be neglected, so that the number of cut-on modes is effectively controlled by R2 ðXÞ and U 0 ðXÞ. As one would expect, the reduction of R2 ðXÞ along a contracting duct has the effect of reducing the right hand side of (45), and therefore tends to cut modes off. However, at the same time as R2 ðXÞ is decreasing, U 0 ðXÞ is increasing through the contraction, and this has the opposite effect of increasing the right hand side of (45), and therefore tends to cut modes on. Of course, for the vortical mean flow used here the precise condition (45) cannot
100 80 60 40 20 0 −20 −40 −60 −80 −100 −35 −30 −25 −20 −15 −10
−5
0
5
10
Fig. 3. The variation of the eigenvalue spectrum in the complex k plane for the expanding duct (13), starting at X ¼0 (open circle symbol) up to X¼ 1.2 (cross symbol). Here we have ω ¼ 26, m¼ 20, U max ¼ 0:5, s ¼ 8, with hard walls.
3944
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
150 100 50 0 −50 −100 −150 −120
−100
−80
−60
−40
−20
0
20
Fig. 4. The variation of the eigenvalue spectrum in the complex k plane for the contracting duct (12), starting at X ¼0 (open circle symbol) up to X ¼1.6 (cross symbol). Other conditions as in Fig. 3.
1
2
0.9 0.8
1.8
0.7 0.6
1.6
0.5
1.4
0.4 0.3
1.2
0.2 0.1
1 0
0.2
0.4
0.6 X
0.8
1
1.2
0
0.2
0.4
0.6
0.8 X
1
1.2
1.4
1.6
Fig. 5. The variation of the square of the amplitude function, jNðXÞ2 j, for the first radial order propagating upstream (solid lines) and downstream (dashed lines) for (a) expanding and (b) contracting ducts. Conditions as in Figs. 3 and 4.
be applied, but the same two physical processes are at work—for the parameter values used in Fig. 4 it is clear that the effect of accelerating mean flow dominates the geometric effect of the narrowing duct, leading to more cut-on modes in the narrower section of the duct. Having determined the axial variation of the eigenmodes, we are now able to solve Eq. (44) to determine the axial variation of the modal amplitude jNðXÞj for each mode, and results are given in Fig. 5 for the first upstream- and downstream-propagating modes. In each case the modal amplitude is set to unity at the nominal starting value of X (so for instance in the contracting case we set Nð1:6Þ ¼ 1 for the upstream-going mode). When given an actual engine geometry with numerical values for the radii R1,2 ðxÞ as functions of the axial coordinate x (recall that x is made non-dimensional by the outer radius at x ¼0), one is faced with the question of what precise value ϵ actually takes, and then whether or not that value is small enough for the asymptotic theory to give a reasonable result. This issue has been considered by Rienstra and Eversman [12], who compare an asymptotic solution assuming slow axial variation with a finite element numerical solution which does not assume slow axial variation, in the case of irrotational mean flow. Good agreement between the two methods is found. They point out that ϵ is a measure of the order of magnitude of the axial rate of change of the duct radius, and that it can therefore be estimated in a given case by evaluating dR1,2 =dx. If we take the ducts described by (12) and (13), and suppose that the axial coordinate x runs from 0 to 1.6 (i.e. the length of the duct is 1.6 times the outer radius at x¼0), then the aspect ratio of the ducts in physical (r,x) space is the same as shown in Fig. 2 in (r,X) space. It then follows that dR1,2 =dx is simply given by directly replacing X by x in (12) and (13) and differentiating with respect to x, so that jdR1,2 =dxj takes values between about 0.07 and 0.1 along the duct. Our small parameter ϵ can therefore be estimated to be 0.1. Since the value ϵ ¼ 0:1 is, if anything, a little smaller than what was found
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3945
in Rienstra and Eversman's case (and they achieved good agreement between their asymptotic and numerical predictions), there is perhaps reason for cautious optimism about the accuracy of the results in this paper (the error involved in the asymptotic solution is formally Oðϵ2 Þ, so only a couple of percent for this value of ϵÞ). However, proper validation against a full numerical code in shear flow is certainly still required. 5. Conclusions In this paper we have described the propagation of acoustic waves along a slowly varying duct carrying radially sheared mean flow. This represents a special case of the analysis presented in [18], where the mean flow had nonzero swirl. By setting the mean swirl to zero, while retaining the more general axial flow profile, we have been able to simplify the calculation of both the mean and unsteady flows compared to [18]. As well as the acoustic modes described in this paper, hydrodynamic modes are also present in the duct. In [19] a method for computing the evolution of these disturbances in a straight duct was presented, using asymptotic analysis in the limit of large azimuthal wavenumber. The extension of that analysis to the question of how the hydrodynamic modes evolve in a slowly varying duct has been considered in [23], but further work is required. Appendix A. Various operators In this appendix we give explicit formulations for the operators introduced in Section 3. The operators appearing in Eqs. (31) and (32) are 0 1 P − 2ωU2 02D0 Dr0 þ D0 ∂r∂ 0 C β 0 0 B C B C 1 B 0 0 0C 2 β B C, 0 L1 ¼ B C B 0 C 0 ω 0 @ A ∂U 0 ∂ ∂r ∂r
where the operator P is defined as P ¼ D0 and
∂U 0 ∂r
0
ω
! ∂2 D0 ∂ ω2 m 2 þ D0 , þ − r ∂r ∂r 2 C 20 r 2 0
0
B B 1 L2 ¼ B B − ∂U 0 @ ∂r 0
D0
0
0
0
0
U0
0
0
The components of the vector f ¼ ðf 1 ,f 2 ,f 3 ,f 4 Þ in Eq. (32) are ∂ i ∂ ðiH 0 D0 Þ þ f1 ¼i ∂X A0 ∂X
U 0 ðω−kU 0 Þ C 20
!
1
C 0 C C: 0 C A
(A.3)
U0 !
D0 A20
þk
D0
(A.2)
i ∂ þ rA0 ∂r
rD0 V 1 ðω−kU 0 ÞA20 C 20
f 2 ¼ 0, f 3 ¼ iV 1
0
,
(A.4) (A.5) (A.6)
∂U 0 ∂H 0 ∂H 0 −V 1 −U 0 : ∂X ∂r ∂X
(A.7)
P
B B 2ωU 0 D0 n B − 2 2 −k D0 C 0 β0 L ¼B B B −D0 ∂r∂ @ †
n
!
∂R0 ∂R0 ∂V 1 ∂A0 ∂U 0 þ iU 0 þ iR0 −i , ∂r ∂X ∂r ∂X ∂r
f 4 ¼ −H 0 The adjoint operator L† is
(A.1)
−k D0
−k
n
k
n ∂U 0
∂r
1 β20
0
0
ω−k U 0
0
0
n
1 ∂U 0 ∂ 0 − 1r ∂r∂ r ∂U ∂r − ∂r ∂r C C C 0 C: C ∂U 0 C A ∂r n ω−k U 0
(A.8)
References [1] D. Pridmore-Brown, Sound propagation in a fluid flowing through an attenuating duct, Journal of Fluid Mechanics 4 (1958) 393–406. [2] W. Eversman, Theoretical models for duct acoustic propagation and radiation, in: H.H. Hubbard (Ed.), Aeroacoustics of Flight Vehicles: Theory and Practice, NASA, 1991, pp. 101–163.
3946
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
[3] P. Doak, P. Vaidya, Attenuation of plane wave and higher order mode sound propagation in lined ducts, Journal of Sound and Vibration 12 (2) (1970) 201–224. [4] G. Vilenski, S.W. Rienstra, Numerical study of acoustic modes in ducted shear flow, Journal of Sound and Vibration 307 (2007) 610–626. [5] G. Vilenski, S.W. Rienstra, On hydrodynamic and acoustic modes in a ducted shear flow with wall lining, Journal of Fluid Mechanics 583 (2007) 45–70. [6] C. Tam, L. Auriault, The wave modes in ducted swirling flows, Journal of Fluid Mechanics 371 (1998) 1–20. [7] V. Golubev, H. Atassi, Sound propagation in an annular duct with mean potential swirling flow, Journal of Sound and Vibration 198 (5) (1996) 601–616. [8] V. Golubev, H. Atassi, Acoustic vorticity waves in swirling flows, Journal of Sound and Vibration 209 (2) (1998) 203–222. [9] C. Heaton, N. Peake, Algebraic and exponential instability of inviscid swirling flow, Journal of Fluid Mechanics 565 (2006) 279–318. [10] N. Peake, A. Parry, Modern challenges facing turbomachinery aeroacoustics, Annual Review of Fluid Mechanics 44 (2012) 227–248. [11] S. Rienstra, Sound transmission in slowly varying circular and annular lined ducts with flow, Journal of Fluid Mechanics 380 (1999) 279–296. [12] S. Rienstra, W. Eversman, A numerical comparison between the multiple-scales and finite-element solution for sound propagation in lined flow ducts, Journal of Fluid Mechanics 437 (2001) 367–384. [13] N. Ovenden, S. Rienstra, Mode matching strategies in slowly varying engine ducts, AIAA Journal 42 (9). [14] E. Brambley, N. Peake, Sound transmission in strongly curved slowly varying cylindrical ducts with flow, Journal of Fluid Mechanics 596 (2008) 387–412. [15] N. Peake, A. Cooper, Acoustic propagation in ducts with slowly varying elliptic cross-section, Journal of Sound and Vibration 243 (3) (2001) 381–401. [16] A. Cooper, Effects of mean entropy on unsteady disturbance propagation in a slowly varying duct with mean swirling flow, Journal of Sound and Vibration 291 (3-5) (2006) 779–801. [17] M. Oppeneer, W. Lazeroms, S. Rienstra, P. Sijtsma, B. Mattheij, Acoustic modes in a duct with slowly varying impedance and non-uniform mean flow and temperature, Proceedings of 17th AIAA/CEAS Aeroacoustics Conference, AIAA-2011-2871. [18] A.J. Cooper, N. Peake, Propagation of unsteady disturbances in a slowly varying duct with mean swirling flow, Journal of Fluid Mechanics 445 (2001) 207–234. [19] A. Cooper, N. Peake, Upstream-radiated rotor–stator interaction noise in mean swirling flow, Journal of Fluid Mechanics 523 (2005) 219–250. [20] A. Cooper, N. Peake, Trapped acoustic modes in aeroengine intakes with swirling flow, Journal of Fluid Mechanics 419 (2000) 151–175. [21] M. Myers, On the acoustic boundary condition in the presence of flow, Journal of Sound and Vibration 71 (3) (1980) 429–434. [22] M. Khorrami, A Chebyshev spectral collocation method using a staggered grid for the stability of cylindrical flows, International Journal for Numerical Methods in Fluids 12 (9) (1991) 825–833. [23] A. Lloyd, Rotor–Stator Interaction Noise, PhD Thesis, University of Cambridge, July 2009. [24] J. Boyd, Chebyshev and Fourier Spectral Methods, Dover, 1999. [25] C. Heaton, Centre modes in inviscid swirling flows and their application to the stability of the Batchelor vortex, Journal of Fluid Mechanics 576 (2007) 325–348.
Contents lists available at SciVerse ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
The propagation of acoustic waves in a slowly varying duct with radially sheared axial mean flow A.E.D. Lloyd, N. Peake n Department of Applied Mathematics & Theoretical Physics, University of Cambridge, Wilberforce Road, Cambridge CB3 0WA, United Kingdom
a r t i c l e i n f o
abstract
Article history: Accepted 1 February 2013 Handling Editor: P. Joseph
We consider the propagation of acoustic waves along a cylindrical duct carrying radially sheared axial mean flow, in which the duct radius is allowed to vary slowly along the axis. In previous work [A.J. Cooper & N. Peake, Journal of Fluid Mechanics 445 (2001) 207–234.] radially sheared axial mean flow with nonzero swirl in a slowly varying duct was considered, but in this paper we set the swirl to zero, thereby allowing simplification of the calculations of both the mean and unsteady flows. In this approach the acoustic wavenumber and corresponding eigenfunction are determined locally, while the wave amplitude is found by solving an evolution equation along the duct. Sample results are presented, including a case in which, perhaps surprisingly, the number of cut-on modes increases as the duct radius decreases. & 2013 Elsevier Ltd. All rights reserved.
1. Introduction There has been considerable interest in the propagation of acoustic waves in flow-bearing ducts for many years, motivated by a number of different industrial applications. In the context of turbomachinery noise, one is often interested in the two different situations; one situation in which the mean flow is both radially sheared and contains a significant swirling component, as found between a rotor row and a stator row; and a second situation in which the swirl is negligible but the axial flow is still radially sheared, as found in the bypass duct downstream of the outlet guide vanes. In both the cases the presence of nonzero mean vorticity is significant and must be included in the modelling. Also in both the cases, but especially in the latter case, the axial variation of the duct cross-section should also be considered. The propagation of acoustic waves in a straight duct carrying axial flow has received much attention, including the seminal work of Pridmore-Brown [1]. Many other early references are given in [2]. A common issue throughout has been the role of wall linings in attenuating the sound, and we would like to mention here in particular an ingenious early paper by Doak and Vaidya [3], in which an approximate analytical expression for attenuation coefficients, valid in a range of geometries, was derived. More recently, Vilenski and Rienstra [4,5] have completed a detailed numerical study over a wide range of frequencies. Although acoustic (pressure-dominated) waves are of primary importance, the sheared axial flow also allows the propagation of hydrodynamic (vorticity-dominated) waves, described in inviscid flow by the critical layer. The inclusion of mean swirl has been considered, for instance by [6–9]; swirl has a significant quantitative effect on the acoustic field, and also gives rise to an additional family of hydrodynamic, ‘nearly convected’ modes located close to the critical layer in wavenumber space. Full details of these hydrodynamic features are given in the above references, and a recent review is provided in [10]. Axially varying ducts have been considered, and in particular the work of Rienstra [11] on slowly varying ducts with radially uniform axial mean flow has proved very influential - see [12–17] for a series of extensions. n
Corresponding author. Tel.: þ 44 1223 339058. E-mail address: [email protected] (N. Peake).
0022-460X/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.jsv.2013.02.038
3938
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
r
r=R2 r=R1
θ X
Fig. 1. Schematic of slowly varying duct.
The particular extension of [11] that we wish to highlight here, however, is the case of acoustic wave propagation in a slowly varying duct containing both mean swirl and radially sheared axial mean flow, which was completed by [18] (with ramifications for noise generation presented in [19]). The analysis in [18] is necessarily complicated, due to the rather general nature of the vortical mean flow considered. The bypass-duct situation of zero swirl but radially sheared axial flow is certainly contained as a special case of the theory in [18], but is not explored there. In fact, by setting the swirl to zero in [18] it turns out that the determination of both the mean and the unsteady flow is simplified, and this is the subject of the current paper. In Section 2 we show how the evolving mean flow can be calculated along the duct, and present the unsteady flow equations; in Section 3 we show how to determine the evolution of the acoustic modes along the slowly varying duct, and in Section 4 we give sample results. In many situations it can be expected that the number of cut-on modes decreases as the duct radius decreases (see for instance [20]), but we present a case in which the number of cut-on modes increases as the duct radius decreases. 2. Problem formulation Consider an annular duct with circular cross section, described using a cylindrical coordinate system ðx,r,θÞ, as shown in Fig. 1. The duct carries a compressible, perfect, isentropic flow, which is made up of a steady mean component plus a small unsteady perturbation. Throughout this paper we non-dimensionalise all lengths by the outer duct radius at x¼0, and velocities and densities by the mean sound speed and mean density on the outer wall at x ¼0, respectively. The duct crosssection varies slowly in the axial direction and we introduce a slow axial scale, X ¼ ϵx, where ϵ⪡1 is a small parameter defined as a measure of the axial slope of the duct walls. The duct inner and outer radii are denoted by R1,2 ðXÞ. The flow has total velocity v, pressure p, density ρ and sound speed c, and is decomposed in the form ~ c~ , ρ, ~ ~ pÞðX,r,θ,tÞ, ðv,c,ρ,pÞðX,r,θ,tÞ ¼ ðV,C,D,PÞðX,rÞ þðv,
(1)
~ c~ , ρ, ~ is the small-amplitude (linear) unsteady perturbation. We will describe how ~ pÞ where ðV,C,D,PÞ is the steady flow and ðv, to determine these two components in the next two subsections. 2.1. Mean flow We consider a mean flow field which has zero swirl but which has axial and radial components which vary with both X and r. It then follows that the mean flow has the asymptotic form ðU,V,C,D,PÞðX,r; ϵÞ ¼ ðU 0 ðX,rÞ,ϵV 1 ðX,rÞ,C 0 ðX,rÞ,D0 ðX,rÞ,P 0 ðX,rÞÞ þ Oðϵ2 Þ,
(2)
while the assumption of a homentropic perfect gas relates P0 and C0 to D0 via P0 ¼
1 γ D γ 0
ðγ−1Þ=2
C 0 ¼ D0
and
:
(3)
At OðϵÞ the continuity equation is ∂ðD0 U 0 Þ 1 ∂ðrD0 V 1 Þ þ ¼ 0, ∂X r ∂r
(4)
and can be satisfied by writing velocity components U0 and V1 in terms of the streamfunction ψðX,rÞ, such that U0 ¼
1 ∂ψ , rD0 ∂r
V1 ¼ −
1 ∂ψ : rD0 ∂X
(5)
Bernoulli's relation (enthalpy is constant on streamlines) can now be written as γ−1
1 2 D0 U þ ¼ HðψÞ, 2 0 γ−1
(6)
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3939
where the enthalpy H is a function of ψ specified by the conditions upstream. For homentropic flow Crocco's relation is
where the mean vorticity ξ is
and taking the r component of (7) gives
V ξ ¼ ∇H,
(7)
∂U 0 ∂V 1 þ ϵ2 eθ , ξ¼ − ∂r ∂X
(8)
∂ 1 ∂ψ ¼ rD0 H′ðψÞ: ∂r rD0 ∂r
(9)
We note here that our non-swirling mean flow Eqs. (6) and (9) take a simpler form than their nonzero-swirl counterparts, Eqs. (17) and (21) in [18]. Now differentiating (6) with respect to r and substituting for H′ðψÞ from (9) we see that ∂D0 =∂r ¼ 0, and therefore that the mean density, and hence from (3) also the mean pressure and sound speed, are functions of X only. The boundaries of the duct are streamlines of the mean flow, along which ψðX,rÞ is constant; the values of these constants will be specified as a part of the initial conditions at X ¼0. For convenience the axial velocity profile at X¼ 0 is written as U 0 ðX ¼ 0,rÞ ¼ g′ðrÞ=r, and the function g(r) is chosen to give a suitable initial profile. Using this profile and integrating Eq. (5) we find ψðX ¼ 0,rÞ ¼ gðrÞ, with the constant of integration taken to be zero without loss of generality (note that our nondimensionalisation implies that D0 ð0Þ ¼ 1). We make the assumption that there is a streamline connecting each point downstream with a point at X¼0 (i.e. no regions of reverse mean flow), allowing us to equate the Bernoulli quantity given by Eq. (6) at X¼0 to the value on the same streamline at an arbitrary distance downstream. After rearrangement this gives rise to the equation ! 1=2 ∂ψ g′ðr s Þ2 2 γ−1 ¼ D0 ðXÞr 1−D0 ðXÞ þ , (10) ∂r γ−1 r 2s where the radial coordinate r s ðX,ψÞ is such that the streamline passing through the point (r,X) also passes through the point ðr s ,0Þ. We now integrate (10) to give " # Z ψ g′ðr Þ2 −1=2 2D20 r 2 −R21 ðXÞ s : (11) þ 1−Dγ−1 dψ ¼ 0 r 2 γ−1 s gðR1 ð0ÞÞ We first set r ¼ R2 ðXÞ in (11) (so that ψ-gðR2 ð0ÞÞ), giving an implicit equation for D0 ðXÞ which can be solved iteratively by marching downstream from X¼ 0. Once D0 ðXÞ has been determined, ψðr,XÞ can be found iteratively from (11), and then finally the mean flow determined from (5). Note that in [18] the calculation of the mean flow with nonzero swirl is more complicated than this: in nonzero swirl the mean pressure, density and sound speed depend on r as well as X; determination of the streamfunction as a function of position requires solution of two coupled ordinary differential equations; and while in principle the method presented in [18] is valid for general upstream conditions, for simplicity only a straightforward case of solid-body rotation and uniform mean flow upstream were considered there. As specific examples we consider two annular ducts, one slowly contracting and the other slowly expanding, based on the examples used in [18]. For the contracting case the duct boundaries are defined by R1 ðXÞ ¼ 0:4482 þ0:05 tanhð2X−2Þ,
R2 ðXÞ ¼ 0:9518−0:05 tanhð2X−2Þ,
(12)
and for the expanding case the duct boundaries are defined by R1 ðXÞ ¼ 0:45180−0:05 tanhð2X−2Þ,
R2 ðXÞ ¼ 1:04820 þ0:05 tanhð2X−2Þ:
(13)
We note that the duct mean line, r ¼ ðR1 ðXÞ þ R2 ðXÞÞ=2, is straight in both the cases. Clearly, other examples are possible in which duct area variation is combined with curvature of the duct mean line, but these additional cases will not be considered here. Following [4] we model the axial velocity profile at X¼ 0 as a parabola, chosen so that the maximum velocity lies in the centre of the duct half way between the inner and outer walls U 0 ðX ¼ 0,rÞ ¼ U max
1 þ sðr−R1 ð0ÞÞð1−rÞ 1 þsð1−R1 ð0ÞÞ2 =4
,
(14)
and where U max ¼ 0:5 and s ¼ 8, corresponding to U 0 ðR1 ð0ÞÞ ¼ 0:76U max . This choice of U0 corresponds to gðrÞ ¼ U max
6r 2 ð1−R1 ð0ÞsÞ þ4sð1 þ R1 ð0ÞÞr 3 −3sr 4 12 þ 3sð1−R1 ð0ÞÞ2
:
(15)
In Fig. 2 we show the mean flow for the contracting duct described in Eq. (12). Note that V1 is small, and therefore that the radial velocity, ϵV 1 , is of negligible magnitude in comparison to the axial component U0. However, the figure also shows that U0 varies significantly along the contracting duct, in which a significant acceleration downstream is accompanied by a reduction in D0.
3940
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
Fig. 2. Leading order components of (a) the axial mean velocity, U 0 ðX,rÞ; (b) the radial mean velocity, V 1 ðX,rÞ. The bold lines denote the duct walls. The axial variation of (c) the mean density, D0 ðXÞ; (d) the mean sound speed C 0 ðXÞ. All for a contracting duct with an initial parabolic mean flow profile (from the left). In (a) the contours are plotted at intervals of 0.05 and in (b) at intervals of 0.005.
2.2. Unsteady flow equations Having obtained the solution for the mean flow we now consider the unsteady disturbance field. Following the work of [7] we decompose the perturbation velocity into a sum of vortical and potential components v~ ¼ uR þ∇ϕ:
(16)
This decomposition is made unique by requiring that the unsteady pressure is expressed solely in terms of the unsteady potential p~ ¼ −D0
Dϕ , Dt
(17)
where D=Dt is the convective derivative D=Dt ¼ ∂=∂t þV⋅∇. The Euler and continuity equations now become DuR þðuR ∇ÞV ¼ −ξ ∇ϕ, Dt ! D 1 D 1 1 − ∇ ðD0 ∇Þ ϕ ¼ ∇ D 0 uR , Dt C 20 Dt D0 D0
(18)
(19)
where ξ is the mean vorticity as defined in (8). In the absence of mean vorticity equation (18) decouples from (19), and uR can be determined everywhere by integrating (19) downstream starting from given conditions at upstream infinity. However, in the present case the unsteady vorticity is fully coupled to the unsteady pressure through the nonzero mean vorticity, and Eqs. (18) and (19) must be solved simultaneously. We will impose boundary conditions to model the presence of an acoustic liner on the inner and outer duct walls of given impedance Z1 and Z2 respectively. For an arbitrary curved wall the correct boundary condition was first given by [21], and has been implemented by [11,18]. It is p~ iω v~ nj ¼ iωþ V ∇−nj ⋅ðnj ⋅∇VÞ , (20) Zj
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3941
where nj , j ¼ 1,2 are the outward (i.e. pointing out of the fluid) unit normal vectors at the wall. For the slowly varying duct these normals are n1 ¼ −
n2 ¼
er −ϵR′1 ex ð1 þϵ2 R′21 Þ1=2 er −ϵR′2 ex
ð1 þ ϵ2 R′22 Þ1=2
,
:
(21)
(22)
Hard-wall boundary conditions are given by the limit Z j -∞. The unsteady flow equations permit two essentially distinct families of solutions, corresponding to acoustic modes which are pressure-dominated and hydrodynamic motions which travel close to the local mean flow speed. In this paper we consider only the acoustic modes. 3. Acoustic modes As in [18] we use the method of multiple scales and seek a solution of the form ðϕ,uRx ,uRr ,uRθ Þðx,r,θ,t; ϵÞ ¼ ðA,H,R,TÞðX,r; ϵÞE:
(23)
Here Z i E ¼ exp iωt−imθ− ϵ
X
kðηÞ dη ,
(24)
ω is the frequency of the acoustic perturbation, m is the azimuthal wavenumber and k(X) is the local axial wavenumber. The amplitudes and the axial wavenumber evolve along the slowly varying duct and therefore depend on the slowly varying X coordinate. Substituting Eq. (23) into Eqs. (18) and (19) gives rise to four coupled equations, which, up to and including OðϵÞ terms, become ! " ! !# ! ∂2 A 1 ∂A Λ2 m2 2 ∂R imT ϵ ∂ iϵ ∂ U0Λ 1 ∂ rD0 V 1 ΛA2 2 , þR þ − þ þ − −k Aþ ðD0 H Þ−ikH ¼ þk D0 A þ r ∂r ∂r r D0 ∂X D0 A ∂X r ∂r ∂r 2 C 20 r 2 C 20 C 20 (25) ∂U 0 ∂R ∂V 1 ∂R ∂A ∂U 0 þR − ¼ −ϵ V 1 þU 0 , iΛRþ ikA ∂r ∂X ∂X ∂r ∂r ∂r
(26)
∂T V 1 T ∂T þ þ U0 , iΛT ¼ −ϵ V 1 ∂r r ∂X
(27)
∂A ∂U 0 ∂U 0 ∂H ∂H þU 0 , ¼ −ϵ H þV1 iΛH þ R þ ∂r ∂r ∂r ∂X ∂X
(28)
where Λ ¼ ω−kU 0 . Similarly, the boundary conditions (20) give ! dRj ∂A D0 Λ2 A iϵ ∂ ∂ ∂V 1 dRj ∂U 0 D0 ΛA2 þR 7 U0 þV1 − ðkA þiHÞ 7 þ iω ¼ ϵω ∂r Zj A ∂X ∂r ∂r dX dX ∂r Zj
(29)
on r ¼ Rj ðXÞ, j¼1,2, where 7 refers to R1 and R2. Eqs. (25)–(29) agree with the corresponding Eqs. (36)–(40) in [18] when the mean swirl is set to zero. We now proceed to expand the amplitudes in (23) in powers of ϵ ðA,H,R,TÞðX,r; ϵÞ ¼ ðA0 ,H 0 ,R0 ,T 0 ÞðX,rÞ þϵðA1 ,H 1 ,R1 ,T 1 ÞðX,rÞ þ Oðϵ2 Þ,
(30)
and substituting into Eqs. (36)–(39) we find at Oð1Þ and OðϵÞ Lψ 0 ¼ 0,
(31)
Lψ 1 ¼ f,
(32)
respectively, where L ¼ L1 −kL2 , and the eigenfunctions take the form ψ n ¼ ðAn ,ηn ,Rn ,iH n Þ for n¼0,1, with ηn ¼ kβ20 An and β20 ¼ 1−U 20 =C 20 . Explicit expressions for L1 , L2 and f are given in Appendix A. We note that the azimuthal component of the perturbation vortical velocity is always equal to zero (i.e. T 0 ¼ T 1 ¼ 0); in the case of nonzero mean swirl [18] all the three components of the vortical velocity are nonzero and must be determined. Eqs. (31) and (32) are to be solved subject to the boundary conditions ∂A0 D0 ðω−kU 0 Þ2 iω þR0 7 A0 ¼ 0 (33) ∂r Zj
3942
and
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
! dRj ∂A1 D0 ðω−kU 0 Þ2 A1 i ∂ ∂ ∂V 1 dRj ∂U 0 D0 ðω−kU 0 ÞA20 þV1 − þ R1 7 ðkA0 þiH 0 Þ 7 þ ¼ω U0 iω A0 ∂X ∂r ∂r ∂r Zj dX dX ∂r Zj
(34)
respectively, which are obtained at OðϵÞ and Oðϵ2 Þ from (29). The leading-order problem defined by (31) and (33) is solved for the local axial wavenumber k(X) as an eigenvalue, with corresponding eigenfunction ψ^ 0 ðX,rÞ. At this order the solution is therefore known up to a slowly varying amplitude function, N(X), so that ψ 0 ðX,rÞ ¼ NðXÞψ^ 0 ðX,rÞ:
(35)
To determine N(X) we will take the inner product of (32) with the homogeneous solution of the adjoint problem. The inner product between the vectors A ¼ ðA1 ,A2 ,A3 ,A4 Þ and B ¼ ðB1 ,B2 ,B3 ,B4 Þ is defined as Z R2 4 ∑ Ann Bn r dr, (36) 〈A,B〉 ¼ R1
n¼1
and we denote the solution of the adjoint problem as ψ †0 ¼ ðY 1 ,Y 2 ,Y 3 ,Y 4 Þ. With this choice of inner product the adjoint problem is found to be L† ψ †0 ¼ 0,
(37)
†
where the adjoint operator L is defined in Appendix A, subject to the boundary condition " # ∂Y n1 D0 ðω−kU 0 Þ2 Y n1 n ∂U 0 Y4 ¼ 0: −D0 7 ∂r ∂r iωZ j
(38)
We therefore have the two systems of equations to solve, one for the eigenvector given in (31) and (33), the other for its adjoint given in (37) and (38). By comparing the two problems, we find that the adjoint eigenvector ψ †0 can be related to the direct eigenvector ψ 0 by ! 2ωU 0 iD0 Hn D0 Rn0 n : (39) þ k An0 , Y 3 ¼ − 0 , Y 4 ¼ − Y 1 ¼ An0 , Y 2 ¼ D0 β20 2 2 ∂U ∂U 0 0 C 0 β0 ∂r ∂r The inner product of (32) with ψ †0 gives 〈ψ †0 ,f〉 ¼ 0, and after some algebra this leads to the governing equation for N(X) in the form FðXÞ
d ðN2 ðXÞÞ ¼ N2 ðXÞGðXÞ: dX
Here F(X) and G(X) are defined as Z FðXÞ ¼ and
R2 ðXÞ
r R1 ðXÞ
U0 Λ C 20
2
! þ k D0 A^ 0
(40)
2
3 2 3 2 2 ^ 02U0 ^ 02U0 D R Λ A D R Λ A 2 1 5 þ4 0 5 dr þ 4 0 iωZ 2 iωZ 1 R2
# n n ∂A^ 0 ^ ^ ^ ^ ^ þ iY 3 f 3 r þ iY 4 f 4 r dr r irD0 H 0 ∂X C 20 R1 ðXÞ R1 ðXÞ 2 0 13 2 0 13 ^2 ^2 D R ∂ ∂ ∂V dR ∂U D Λ A D R ∂ ∂ ∂V dR ∂U D Λ A 0 2 1 2 0 0 0 1 1 1 0 0 0 0 @ A 5 −4 @ A5 þV1 − þV1 − U0 þ U0 þ −4 ∂X ∂r ∂r ∂X ∂r ∂r iω dX ∂r Z2 iω dX ∂r Z1 d GðXÞ ¼ − dX
Z
R2 ðXÞ
"
U0 Λ
!
(41) R1
# Z 2 ^ ^ ^ þ k D0 A 0 þ irD0 A 0 H 0 dr þ
R2
R2 ðXÞ
"
(42) R2
respectively. The functions f^ 3 and f^ 4 refer to the components of the vector f given in Appendix A, but with non-hatted variables replaced with their hatted equivalents. The general solution to (40) is now simply Z X GðηÞ dη , (43) N2 ðXÞ ¼ N 2 ð0Þ exp 0 FðηÞ where we have chosen to fix the amplitude relative to its value at X¼0. Eq. (43) provides the most convenient form for evaluation of the variation of N(X). However, for ease of comparison with the radially uniform mean flow solution in [11], we also note that we can rearrange (40) into the form #
Z R2 ðXÞ Z R2 ðXÞ " n n d ∂A^ 0 2 2 d 2 ^ ^ ^ ^ ^ ^ ^ ðFN Þ ¼ −N þ iY 3 f 3 r þ iY 4 f 4 r dr: iD0 A 0 H 0 dr þN irD0 H 0 (44) dX dX ∂X R1 ðXÞ R1 ðXÞ In the special case of potential disturbances in radially uniform mean flow, U 0 ðX,rÞ is replaced by U 0 ðXÞ and the vortical component of the unsteady velocity becomes zero (when the mean vorticity is zero, zero unsteady vorticity upstream
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3943
implies zero unsteady vorticity throughout the flow). The right hand side of Eq. (44) then becomes exactly zero, and (44) agrees exactly with the corresponding equation (4.10) in [11]. 4. Results Our first step is to solve numerically the leading-order eigenvalue problem L1 ¼ kL2 , together with the boundary condition (33), to obtain the unknown axial wavenumbers k and the corresponding eigenfunctions, at each X. As in [18], to do this we employ a Chebyshev spectral collocation method of the form given by [22], in which we choose nr collocation points in the radial direction. Full details of this method can be found in these references and in [23]. One of the well-known characteristics of pseudo-spectral methods is the appearance of spurious eigenvalues, for example see [24]. In our problem these include spurious modes located in an egg-shaped region around the critical layer, as noticed by [25]. Following the work of [14], we apply two filtering procedures to remove these spurious eigenvalues; the first removes eigenvalues whose position varies significantly as nr is varied, while the second removes eigenvalues for which the corresponding eigenfunction is not sufficiently well resolved. Full details are again given in [23]. We also note that if one were interested in identifying hydrodynamic modes, as opposed to the acoustic modes considered throughout this paper, then even more care would have to be taken when separating genuine and spurious eigenvalues. For instance, in work on swirling ducted flow [9], classical numerical shooting was found to be effective, even for modes close to the critical layer. This technique could be included as a further check on the results if the above pseudospectral methods were applied to the hydrodynamic modes. We now look in detail at the axial variation of the spectra for the expanding duct (Fig. 3) and the contracting duct (Fig. 4). In both the cases, there are two cut-on modes propagating downstream and two cut-on modes propagating upstream at X¼0. As we move downstream the cut-on and cut-off eigenvalues move in the complex k plane. In the expanding case there are still just two upstream and two downstream propagating cut-on modes at X ¼1.2, and the modes which were cut-off at X¼0 remain cut-off at X¼ 1.2—in fact this continues to be the case as X-∞. However, more complicated behaviour is apparent in the contracting case, in which the less cut-off modes at X¼0 move significantly towards the real axis downstream, and indeed by X¼1.6 (Fig. 4) two pairs of original cut-off modes have moved onto the real axis, so that there are now four cut-on modes propagating downstream and four cut-on modes propagating upstream. This is perhaps surprising, since one normally thinks of modes becoming cut-off as they propagate along a narrowing duct; see for instance [20], in which a model for acoustic resonance in an intake duct involves trapping between the duct narrowing upstream and the swirl behind the fan downstream. To understand this behaviour, we consider the special case in which the duct is hollow (i.e. R1 ðXÞ ¼ 0) and the flow far upstream is uniform (i.e. s ¼ 0 and U 0 ðX,rÞ ¼ U 0 ðXÞ, as studied by [11]). In this case there are well-known analytic expressions for the axial eigenvalues (see for example [11, Eq. 4.9]), and if we simplify further to the case of hard walls we find that the number of cut-on modes corresponds to the number of zeros of J′m ðzÞ, j′mn for n ¼ 1,2,…, which satisfy j′mn o
ωR2 C 0 ð1−U 20 =C 20 Þ
:
(45)
In (45) the X variation of C 0 ðXÞ is rather weak (see Fig. 2b) and can be neglected, so that the number of cut-on modes is effectively controlled by R2 ðXÞ and U 0 ðXÞ. As one would expect, the reduction of R2 ðXÞ along a contracting duct has the effect of reducing the right hand side of (45), and therefore tends to cut modes off. However, at the same time as R2 ðXÞ is decreasing, U 0 ðXÞ is increasing through the contraction, and this has the opposite effect of increasing the right hand side of (45), and therefore tends to cut modes on. Of course, for the vortical mean flow used here the precise condition (45) cannot
100 80 60 40 20 0 −20 −40 −60 −80 −100 −35 −30 −25 −20 −15 −10
−5
0
5
10
Fig. 3. The variation of the eigenvalue spectrum in the complex k plane for the expanding duct (13), starting at X ¼0 (open circle symbol) up to X¼ 1.2 (cross symbol). Here we have ω ¼ 26, m¼ 20, U max ¼ 0:5, s ¼ 8, with hard walls.
3944
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
150 100 50 0 −50 −100 −150 −120
−100
−80
−60
−40
−20
0
20
Fig. 4. The variation of the eigenvalue spectrum in the complex k plane for the contracting duct (12), starting at X ¼0 (open circle symbol) up to X ¼1.6 (cross symbol). Other conditions as in Fig. 3.
1
2
0.9 0.8
1.8
0.7 0.6
1.6
0.5
1.4
0.4 0.3
1.2
0.2 0.1
1 0
0.2
0.4
0.6 X
0.8
1
1.2
0
0.2
0.4
0.6
0.8 X
1
1.2
1.4
1.6
Fig. 5. The variation of the square of the amplitude function, jNðXÞ2 j, for the first radial order propagating upstream (solid lines) and downstream (dashed lines) for (a) expanding and (b) contracting ducts. Conditions as in Figs. 3 and 4.
be applied, but the same two physical processes are at work—for the parameter values used in Fig. 4 it is clear that the effect of accelerating mean flow dominates the geometric effect of the narrowing duct, leading to more cut-on modes in the narrower section of the duct. Having determined the axial variation of the eigenmodes, we are now able to solve Eq. (44) to determine the axial variation of the modal amplitude jNðXÞj for each mode, and results are given in Fig. 5 for the first upstream- and downstream-propagating modes. In each case the modal amplitude is set to unity at the nominal starting value of X (so for instance in the contracting case we set Nð1:6Þ ¼ 1 for the upstream-going mode). When given an actual engine geometry with numerical values for the radii R1,2 ðxÞ as functions of the axial coordinate x (recall that x is made non-dimensional by the outer radius at x ¼0), one is faced with the question of what precise value ϵ actually takes, and then whether or not that value is small enough for the asymptotic theory to give a reasonable result. This issue has been considered by Rienstra and Eversman [12], who compare an asymptotic solution assuming slow axial variation with a finite element numerical solution which does not assume slow axial variation, in the case of irrotational mean flow. Good agreement between the two methods is found. They point out that ϵ is a measure of the order of magnitude of the axial rate of change of the duct radius, and that it can therefore be estimated in a given case by evaluating dR1,2 =dx. If we take the ducts described by (12) and (13), and suppose that the axial coordinate x runs from 0 to 1.6 (i.e. the length of the duct is 1.6 times the outer radius at x¼0), then the aspect ratio of the ducts in physical (r,x) space is the same as shown in Fig. 2 in (r,X) space. It then follows that dR1,2 =dx is simply given by directly replacing X by x in (12) and (13) and differentiating with respect to x, so that jdR1,2 =dxj takes values between about 0.07 and 0.1 along the duct. Our small parameter ϵ can therefore be estimated to be 0.1. Since the value ϵ ¼ 0:1 is, if anything, a little smaller than what was found
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
3945
in Rienstra and Eversman's case (and they achieved good agreement between their asymptotic and numerical predictions), there is perhaps reason for cautious optimism about the accuracy of the results in this paper (the error involved in the asymptotic solution is formally Oðϵ2 Þ, so only a couple of percent for this value of ϵÞ). However, proper validation against a full numerical code in shear flow is certainly still required. 5. Conclusions In this paper we have described the propagation of acoustic waves along a slowly varying duct carrying radially sheared mean flow. This represents a special case of the analysis presented in [18], where the mean flow had nonzero swirl. By setting the mean swirl to zero, while retaining the more general axial flow profile, we have been able to simplify the calculation of both the mean and unsteady flows compared to [18]. As well as the acoustic modes described in this paper, hydrodynamic modes are also present in the duct. In [19] a method for computing the evolution of these disturbances in a straight duct was presented, using asymptotic analysis in the limit of large azimuthal wavenumber. The extension of that analysis to the question of how the hydrodynamic modes evolve in a slowly varying duct has been considered in [23], but further work is required. Appendix A. Various operators In this appendix we give explicit formulations for the operators introduced in Section 3. The operators appearing in Eqs. (31) and (32) are 0 1 P − 2ωU2 02D0 Dr0 þ D0 ∂r∂ 0 C β 0 0 B C B C 1 B 0 0 0C 2 β B C, 0 L1 ¼ B C B 0 C 0 ω 0 @ A ∂U 0 ∂ ∂r ∂r
where the operator P is defined as P ¼ D0 and
∂U 0 ∂r
0
ω
! ∂2 D0 ∂ ω2 m 2 þ D0 , þ − r ∂r ∂r 2 C 20 r 2 0
0
B B 1 L2 ¼ B B − ∂U 0 @ ∂r 0
D0
0
0
0
0
U0
0
0
The components of the vector f ¼ ðf 1 ,f 2 ,f 3 ,f 4 Þ in Eq. (32) are ∂ i ∂ ðiH 0 D0 Þ þ f1 ¼i ∂X A0 ∂X
U 0 ðω−kU 0 Þ C 20
!
1
C 0 C C: 0 C A
(A.3)
U0 !
D0 A20
þk
D0
(A.2)
i ∂ þ rA0 ∂r
rD0 V 1 ðω−kU 0 ÞA20 C 20
f 2 ¼ 0, f 3 ¼ iV 1
0
,
(A.4) (A.5) (A.6)
∂U 0 ∂H 0 ∂H 0 −V 1 −U 0 : ∂X ∂r ∂X
(A.7)
P
B B 2ωU 0 D0 n B − 2 2 −k D0 C 0 β0 L ¼B B B −D0 ∂r∂ @ †
n
!
∂R0 ∂R0 ∂V 1 ∂A0 ∂U 0 þ iU 0 þ iR0 −i , ∂r ∂X ∂r ∂X ∂r
f 4 ¼ −H 0 The adjoint operator L† is
(A.1)
−k D0
−k
n
k
n ∂U 0
∂r
1 β20
0
0
ω−k U 0
0
0
n
1 ∂U 0 ∂ 0 − 1r ∂r∂ r ∂U ∂r − ∂r ∂r C C C 0 C: C ∂U 0 C A ∂r n ω−k U 0
(A.8)
References [1] D. Pridmore-Brown, Sound propagation in a fluid flowing through an attenuating duct, Journal of Fluid Mechanics 4 (1958) 393–406. [2] W. Eversman, Theoretical models for duct acoustic propagation and radiation, in: H.H. Hubbard (Ed.), Aeroacoustics of Flight Vehicles: Theory and Practice, NASA, 1991, pp. 101–163.
3946
A.E.D. Lloyd, N. Peake / Journal of Sound and Vibration 332 (2013) 3937–3946
[3] P. Doak, P. Vaidya, Attenuation of plane wave and higher order mode sound propagation in lined ducts, Journal of Sound and Vibration 12 (2) (1970) 201–224. [4] G. Vilenski, S.W. Rienstra, Numerical study of acoustic modes in ducted shear flow, Journal of Sound and Vibration 307 (2007) 610–626. [5] G. Vilenski, S.W. Rienstra, On hydrodynamic and acoustic modes in a ducted shear flow with wall lining, Journal of Fluid Mechanics 583 (2007) 45–70. [6] C. Tam, L. Auriault, The wave modes in ducted swirling flows, Journal of Fluid Mechanics 371 (1998) 1–20. [7] V. Golubev, H. Atassi, Sound propagation in an annular duct with mean potential swirling flow, Journal of Sound and Vibration 198 (5) (1996) 601–616. [8] V. Golubev, H. Atassi, Acoustic vorticity waves in swirling flows, Journal of Sound and Vibration 209 (2) (1998) 203–222. [9] C. Heaton, N. Peake, Algebraic and exponential instability of inviscid swirling flow, Journal of Fluid Mechanics 565 (2006) 279–318. [10] N. Peake, A. Parry, Modern challenges facing turbomachinery aeroacoustics, Annual Review of Fluid Mechanics 44 (2012) 227–248. [11] S. Rienstra, Sound transmission in slowly varying circular and annular lined ducts with flow, Journal of Fluid Mechanics 380 (1999) 279–296. [12] S. Rienstra, W. Eversman, A numerical comparison between the multiple-scales and finite-element solution for sound propagation in lined flow ducts, Journal of Fluid Mechanics 437 (2001) 367–384. [13] N. Ovenden, S. Rienstra, Mode matching strategies in slowly varying engine ducts, AIAA Journal 42 (9). [14] E. Brambley, N. Peake, Sound transmission in strongly curved slowly varying cylindrical ducts with flow, Journal of Fluid Mechanics 596 (2008) 387–412. [15] N. Peake, A. Cooper, Acoustic propagation in ducts with slowly varying elliptic cross-section, Journal of Sound and Vibration 243 (3) (2001) 381–401. [16] A. Cooper, Effects of mean entropy on unsteady disturbance propagation in a slowly varying duct with mean swirling flow, Journal of Sound and Vibration 291 (3-5) (2006) 779–801. [17] M. Oppeneer, W. Lazeroms, S. Rienstra, P. Sijtsma, B. Mattheij, Acoustic modes in a duct with slowly varying impedance and non-uniform mean flow and temperature, Proceedings of 17th AIAA/CEAS Aeroacoustics Conference, AIAA-2011-2871. [18] A.J. Cooper, N. Peake, Propagation of unsteady disturbances in a slowly varying duct with mean swirling flow, Journal of Fluid Mechanics 445 (2001) 207–234. [19] A. Cooper, N. Peake, Upstream-radiated rotor–stator interaction noise in mean swirling flow, Journal of Fluid Mechanics 523 (2005) 219–250. [20] A. Cooper, N. Peake, Trapped acoustic modes in aeroengine intakes with swirling flow, Journal of Fluid Mechanics 419 (2000) 151–175. [21] M. Myers, On the acoustic boundary condition in the presence of flow, Journal of Sound and Vibration 71 (3) (1980) 429–434. [22] M. Khorrami, A Chebyshev spectral collocation method using a staggered grid for the stability of cylindrical flows, International Journal for Numerical Methods in Fluids 12 (9) (1991) 825–833. [23] A. Lloyd, Rotor–Stator Interaction Noise, PhD Thesis, University of Cambridge, July 2009. [24] J. Boyd, Chebyshev and Fourier Spectral Methods, Dover, 1999. [25] C. Heaton, Centre modes in inviscid swirling flows and their application to the stability of the Batchelor vortex, Journal of Fluid Mechanics 576 (2007) 325–348.