Solving Galbrun’s Equation with a Discontinuous galerkin Finite Element Method

Over many years, scientists and engineers have developed a broad variety of mathematical formulations to investigate the propagation and interactions with flow of flow-induced noise in early-stage of product design and development. Beside established theories such as the linearized
 Euler equations (LEE), the linearized Navier–Stokes equations (LNSE) and the acoustic perturbation equations (APE) which are described in an Eulerian framework, Galbrun utilized a mixed Lagrange–Eulerian framework to reduce the number of unknowns by representing perturbations by
 means of particle displacement only. Despite the advantages of fewer degrees of freedom and the reduced effort to solve the system equations, a computational approach using standard continuous finite element methods (FEM) suff ers from instabilities called spurious modes that pollute the solution.
 In this work, the authors employ a discontinuous Galerkin approach to overcome the difficulties related to spurious modes while solving Galbrun's equation in a mixed and pure displacement based formulation. The results achieved with the proposed approach are compared with results from previous
 attempts to solve Galbrun's equation. The numerical determination of acoustic modes and the identification of vortical modes is discussed. Furthermore, case studies for a lined-duct and an annulus supporting a rotating shear-flow are investigated.


Introduction
Acoustic noise reduction, which is a wide matter of concern in industry, calls for a better understanding of the complex phenomena that occur when an acoustic wave propagates in a mean flow. A majority of research performed in aeroacoustics and computer aided engineering (CAE) has been aimed at the aircraft noise community, arguably starting by the work of Lighthill, [25]. Commonly, aircraft noise research traditionally focusses on high Mach number and high Reynolds number free-field jet flows. In high-speed jets, noise generation is considered to be of quadrupole type, caused by unsteady nonlinear mechanisms. The methodologies developed in CAE have rejected this, in their focus on time-domain solutions of the non-linear Navier-Stokes equations, either as Direct Numerical Simulations (DNS) where no turbulence models are included, to turbulence models such as Reynolds-Averaged Navier-Stokes (RANS) and Large Eddy Simulation (LES) codes. In wall-bound and internal flows at low Mach numbers, the sound generating mechanisms are however governed by fundamentally different physics than that of free-field jet noise. When an air flow is obstructed by a change of geometry, such as a sharp corner or a bifurcation, flow instabilities and vortices are generated. As Received 14 May 2019, accepted 23 September 2019. these vortices impinge on boundaries, sound impulses are generated.
A less explored field of aeroacoustics is that of pure wave propagation in inhomogeneous media with arbitrary mean flows, as this is disconnected from the noise generation processes. The conceptual difference in the simulation of sound generation and sound propagation is large enough to justify a treatment of the two as separate topics. In regions outside of acoustic sources, the acoustic quantities are often small in comparison to the flow-field quantities. In many cases, it can be assumed that the flow-field affects the sound waves, whereas the sound waves do not induce the flow-field. Thus, the perturbations about the mean flow are often small enough to justify linearization. This enables a two-stage treatment of the acoustic wave propagation: firstly, the mean flow can be calculated without the need to consider any acoustic waves, and secondly, the sound waves can be calculated as perturbations about the mean flow-field. Also, as a consequence of the linearization, a frequency domain approach can be taken. A main benefit of a frequency-domain approach, as opposed to a time-domain approach, is the significant reduction of computational time in case of harmonic excitation. Since most research efforts have been aimed at jet noise generation, where unsteady simulations are needed, only a few studies on frequency-domain aeroacoustics are available.
Galbrun [18] proposed a displacement based description in a Lagrange-Eulerian mixed frame for analyzing the propagation of sound waves in moving fluids. Since only the displacement field is the unknown quantity, Galbrun's equation represents a potent alternative to the methods mentioned above. Among the possibility of reducing the degrees of freedom, boundary conditions can naturally be expressed, i.e. in terms of boundary displacement.
Despite the positive aspects of Galbrun's equation, socalled spurious modes exist when extracting the eigenvalues of the associated boundary value problem utilizing a standard, unstabilized finite element method. The same holds for the LNSE and LEE formulation. These spurious modes possibly pollute the solution when using the standard continuous Galerkin discretization [20,13]. Various attempts have been published in the literature for handling the problems associated with spurious modes, cf. [5,35,41,42]. The work presented in this paper aims as a step to the development of a simulation methodology for alternative strategies for flow acoustics.
Bécache et al. [2] and Bonnet et al. [6,5] have demonstrated that the direct displacement-based formulation associated with Galbrun's equation may produce erroneous or spurious solutions if the finite element method is based on simple continuous finite elements. In particular, it is proposed in [5] that a regularized reformulation of the variational equation for uniform and shear flows can produce robust solutions by damping them out. Other authors, such as Dietzsch et al. [12], have used finite element functions which can shift, but not remove, the location of spurious solutions to higher frequencies and higher damping values. In this paper, a formulation based on the displacement variables is presented which minimizes the appearance of spurious solutions without the need of a regularized reformulation.
The paper is organized as follows: Section 2 sets up the Galbrun equations from basic principles which leads to Section 3 a description of the numerical scheme. The authors utilize a discontinuous Galerkin (DG) method, cf. [9,10,1,39], for discretizing Galbrun's equation for which a time-harmonic behavior is assumed. To account for inter-element fluxes, a local Lax-Friedrichs flux, cf. [21], is used and discussed. In order to highlight the benefits of the proposed method, solutions are compared to the state of the art methods such as LNSE and LEE as case studies in Section 3 which includes a discussion on the influence of an appropriate flux factor. Further, a filtering (Lagrange multiplier) technique is applied to exclude non relevant modes from the solution space. Finally, in Section 4 examples including a lined duct and a circulating flow within an annulus are investigated.

Theory and Numerical Method
In this section, the fundamental principles for deriving Galbrun's equation are outlined. For additional insight, the reader is referred to the literature [18,41,33,6,12,19,20,23,34]. In order to give a comprehensive description, some mathematical fundamentals are required. Hereafter, a Cartesian coordinate system defined by the orthonormal directions '1' and '2' is used, such that a two-dimensional space is considered. Further, a vector component description together with Einstein's summation convention is assumed to indicate component summation for repeated indices. Hereafter, the domain of interest is Ω F ⊂ R 2 which is bounded by a closed surface, Γ F .

Lagrange-Eulerian frame
For deriving Galbrun's equation, one has to be familiar with the concept of a mixed Lagrange-Eulerian frame in the view of continuum mechanics and the associated appropriate balance equations. Despite the fact that the theory to derive Galbrun's equation is well published in the literature, the authors present the basic principles for better readability.
Two distinguished states are considered. In the first state namely the perturbed state, any given particle is defined by its spatial coordinate position y(t) where in the second state or the unperturbed state or reference state, the same particle takes the position x(t), cf. Figure 1. The vector components of the Lagrangian displacement w(t) are defined as the difference of these two states, i.e.
Further, any given field quantity Φ takes the form in the Eulerian frame with the Eulerian perturbation Φ (y l , t) and in the Lagrangian frame with the corresponding Lagrangian perturbationΦ(x l , t), cf. Poirée [33]. Since the application of Galbrun's equation is basically a linear perturbation procedure, the mean flow quantities, Φ 0 , can be calculated from an adequate preceding stationary boundary value problem. It must be noted that these mean flow quantities Φ 0 are described specifically in Eulerian coordinates, cf. Brazier [7]. Combining equations (1)-(3) while assuming that the mean values Φ 0 are slowly varying, i.e. associated gradients in time and space are integrable continuous functions. Together with a Taylor expansion up to the linear terms, one identifies the relation between Eulerian and Lagrangian perturbations where (), j = ∇() = ∂() ∂x j represents the Nabla-Operator with the corresponding spatial derivative in the j-direction. From equation (4), it is apparent that if the spatial gradient of the mean values vanishes, the Eulerian and the Lagrangian perturbations are equal, see [31].

Galbrun's equation formulation
The acoustic radiation of a source produces a small perturbation of the physical quantities such as pressure and density. The propagation of that small perturbation is governed by the Galbrun equation which is a linear equation whose unknown w is the perturbation of displacement. The well known conservation equations of fluid dynamics, namely the mass, momentum and energy balance equations in an Eulerian frame are used to derive Galbrun's equation. Assuming small perturbations, the displacement based expression is formulated. Further, it is assumed that the fluid of interest is a perfect inviscid gas with adiabatic thermodynamic properties (i.e. isentropic material behavior). Under these assumptions, Galbrun's equation is stated, cf. the references [18,41,7,31].
In addition, appropriate initial boundary for w j (t = 0) and D(w j )/Dt(t = 0) in Ω F and on Γ F must be stated that fulfill the boundary conditions on Γ F . Further, where ρ 0 represents the mean flow mass density, p 0 the mean flow pressure, c 0 the speed of sound and v 0k the mean flow velocity in the k-direction, respectively. The surface normal vector n on Γ F is pointing outward of the domain Ω F .
Following the arguments considered in the works by Gabard et al. [17], Treyssède et. al [41,40] and Wang et al. [42], a mixed formulation can be achieved by introducing a Lagrangian pressure perturbatioñ Following this, combining equations (5) and equation (8) yields To revert the displacement field back to the well known acoustic pressure p in the associated Eulerian frame, Expressions (4) and (8) are required. Note that the assumptions of slowly varying mean flow quantities such as ρ 0 , p 0 and v 0k must still hold. Further, a time harmonic dependency is assumed to convert the Galbrun equation from the time to the frequency domain, i.e. any given quantity takes the form φ(x, t) = φ (x)e −jωt , while j depicts the imaginary unit, ω = 2πf represents the angular frequency andφ(x) the complex amplitude. Hereafter, the hat symbol is omitted to improve readability. It is thought that all quantities and fields are in a time harmonic regime. Furthermore, the operator of the material time derivative reduces to In equation (6), rigid boundary conditions representing hard walls have been introduced. To account for more complicated boundary conditions, we present an admittance boundary condition for the displacement based Galbrun equation, see [30,29], expressed in Eulerian quantities, where v f is the normal component of the Eulerian fluid velocity perturbation, v s the normal component of the Eulerian structure velocity perturbation,Ȳ = Y/(ρ 0 c 0 ) the boundary admittance and its normalized partȲ and p the Eulerian acoustic pressure perturbation. Applying these and assuming zero structural velocity v s = 0, the equivalent admittance boundary condition for Galbrun's equation may be derived. Making use of equation (4), this leaves a relation between admittance and displacement = Y − c 2 0 ρ 0 w l,l − w l p 0,l . Note thatṽ k = D(w k )/Dt has been used for the relation between the Lagrangian velocity perturbationṽ k , when applying equation (4), and the Lagrangian displacement. If the mean pressure p 0 is constant and the flow velocity v 0k is constant and orthogonal to the surface normal vector n l , equation (12) reduces to a familiar form for boundary admittance which is in agreement with Dietzsch et al. [12].

Numerical method
Discontinuous Galerkin finite element methods (DG-FEM) combine favourable features of finite element methods (FEM) and finite volume methods (FVM) with strong mathematical foundation. DG-FEM possess a number of favourable properties especially in hydrodynamic, uniform and non-uniform flow problems. They are robust and high-order accurate and are able to capture physical phenomena common to mixed finite element problems other methods cannot reach. For this problem we make use of DG-FEM to overcome an inherent issue. Generally, DG-FEM combines the flexibility of introducing high-order FEM schemes with the flexibility of FVM to formulate numerical schemes locally, which can reflect flow, for example. However, finite volume methods are often too inaccurate and diffusive when applied, especially, to wave propagation problems. The basic FVM is a form of the lowest-order DG-FEM. So, it makes sense to increase accuracy but keep the conservation part by developing DG-FEM further.
Considering Figure 2, it is clear that the global approximation u(x) is not well-defined at the boundaries of each element. As Figure 2 demonstrates, at each node x, two solutions exist belonging to the respective adjacent elements so that both, a left element D k and a right element D k+1 evaluate the approximation u(x). Since we do not enforce continuity over the boundaries of the elements as in standard finite element methods, we are not guaranteed that u il (x) = u ir (x), and thus, it appears u(x) is not uniquely defined. To establish a connection between elements, we introduce a flux condition. For this study, a Lax-Friedrich scheme, which blends between central and upwind flux, is chosen. This defines the DG-FEM formulation. For a detailed description of the numerical scheme and possible stabilization techniques, the authors refer to the literature [1,10,9,39].
In this section, the finite element method utilized to solve Galbrun's equation is described. To do this, we first construct the weak form of Galbrun's equation with appropriate test functions that exist in a suitable mathematical space. We define the inner product as in the space of functions denoted by L 2 (Ω F ) where all functions are square integrable over the domain Ω F . In the weak formulation, this inner product is used to define a local inner product and norm such that with where Ω h represents the approximated domain due to the discretization procedure. Figure 2. A Lax-Friedrichs flux is defined across a shared element boundary, shown as a red line Γ between the two elements D k and D k+1 . Displacement across this boundary can be discontinuous, i.e. u il = u ir .
This allows the local unknowns to be discontinuous from one element D k to the other.
To account for inter element fluxes, the Lax-Friedrichs-Flux, cf. Equation (17), is chosen, where the scaling factor α needs to be defined. For Galbrun's equation, the flux term scales with c 2 0 . Therefore, in accordance to Hesthaven and Warburton [21], the flux constant α is set to α = 10 6 , see Section 3.4. Figure 2 illustrates the definition of the Lax-Friedrichsflux between two elements D k and D k+1 , respectively. The flux in element D k across the boundary Γ is denoted by f (u il ) and vice versa from element D k+1 across Γ with f (u ir ) where at the boundary the unknowns can take the value u il in element D k and u ir in element D k+1 . The corresponding outward normal vectors are n k and n k+1 .
Applying all the forgoing principles, the weak form of Galbrun's equation for each element domain Ω e h with boundaries Γ e h reads as with the piecewise discontinuous complex conjugated test functionsw * k . Expanding the material time derivative using equation (7) and simplifying ∂(v 0k )/∂t = 0, ρ 0 ,l = 0 and After integrating by parts and utilizing Green's identity, the weak form can be rearranged as In order to apply the discontinuous Galerkin method, the boundary integral in equation (20) does not vanish and accounts for inter-elemental fluxes. At this stage, any constraint condition such as the restriction to the rotational field of w, e.g. ∇ × w = 0, can be easily integrated using a Lagrange multiplier λ h in the form of where e kmn is known as the permutation tensor. For the remainder of this paper the Lagrange multiplier is discretized using piecewise discontinuous linear basis functions.
To demonstrate how the inter-elemental fluxes are set up, as an example, the last part of the boundary integral in equation (20) is reformulated according to the notation highlighted in Figure 2.
Utilizing the well known approach of partial integration and discretizing the unknowns as well as the test functions according to Galerkin, a matrix formulation in the form of can be achieved. M represents the mass matrix and D , C the damping and stiffness matrices, respectively. The system of equations is analyzed in terms of eigenvalues and eigenvectors using a Krylov subspace method where the size of the subspace is set to be twice as large as the number of eigenvalues of interest, cf. [8,14,37]. For this purpose, the software tool Comsol Multiphysics together with Matlab is used to set up the model as well as the system of equations and finally solve for eigenvalues and eigenvectors.

Finite element model
In this section, we present a parameter study to investigate the performance of applying a discontinuous Galerkin The no-flow case for Galbrun's equation was first studied in the 1970s [20] and was proved to exhibit spurious circulation or rotational modes. In the case spurious modes were solutions of the source-free problem associated with non-zero real-valued frequencies. Determined numerically, these modes could be shifted to higher frequencies by refining the mesh or by enforcing the numerical method to avoid these solutions [12]. It was not until Wang and Bathe [42] showed a mixed pressuredisplacement formulation, using a certain type of finite element, could lead to numerical solutions without spurious modes.
For the no-flow case, a pure displacement formulation can make use of edge elements which are also known as Raviart-Thomas elements [3]. It was shown in [36] that they are of similar efficiency as Lagrangian elements when used for the pure pressure formulation. However, with mean-flow, a mixed formulation of Galbrun's equation can develop numerical spurious modes due to the inherent mathematical nature of the problem [17].

Finite duct with uniform mean flow
Since this article focuses on the finite element method used to compute a solution of the Galbrun equation, we consider the artificial case set in a bounded domain, since analytical results are available. Figure 3 illustrates a rectangular duct configuration in two dimensions with height H = 0.5 m and length L = 3.4 m. The air filled duct is bounded by acoustically hard walls, i.e. the surface normal particle velocity is set to zero. In terms of particle displacement within Galbrun's equation, the following expression is equivalent to the hard walled boundary condition: For this test case, the Mach number is defined by the ratio of the homogeneous constant mean flow velocity v 01 and the constant speed of sound c 0 = 340 m/s. The flow velocity is varied so that the Mach number takes values from 0 to 0.3.  Figure 4 illustrates the finite element mesh consisting of 930 triangular elements. Table I lists all investigated model configurations for the parameter analysis. As an example, the abbreviation "LNSE ρ1v2" is understood as the linearized Navier-Stokes equation where a conventional continuous Galerkin discretization method is used and the basis functions are of 1 st -order (or linear) for the density ρ and of 2 nd -order (or quadratic) for the velocity v. The temperature field within the LNSE and the pressure field within the LEE are omitted since all processes are considered adiabatic. Further, no stabilization scheme is applied. The mixed formulations are solved numerically using Taylor-Hood elements, see [12], which are commonly adopted as basis functions for mixed-formulations.
In order to compare the results of the eigenvalue extraction, analytical solutions of the convected Helmholtz equation are taken as references. For the duct case with plane wave propagation, the eigenfrequencies of the acoustic modes are calculated as follows (cf. Dietzsch et al. [12]).

Convergence study
For a convergence study the duct geometry from the previous numerical model, Section 3.1, was adopted to investigate the numerical scheme. On applying a discontinuous Galerkin method to the displacement based form of Galbrun's equation, Figure 5 illustrates two structured mesh configurations, i.e. a first and a second level of mesh refinement. In total seven mesh refinements, from eight triangular meshes, have been conducted. The relative error for each mesh refinement can be computed since the exact solution is easily derived. Relative errors, ε are measured in discrete L 2 norms and h := min(H x , h y ) > 0 is the size of a triangle element length. To illustrate how the numerical schemes improve with mesh-refinement, as h → 0, experimental convergence rates are illustrated in Figure 6. It is clear that solutions follow an asymptotic convergence rate which suggests numerical stability. The outliers in Figure 6 for the first mesh refinements are thought to be a result of a spatial undersampling of the associated eigenvectors and can thus be ignored. To the authors' knowledge, it is not clear why the convergence rate remains linear when using quadratic elements. Nevertheless due to computational costs, the authors use linear and quadratic elements hereafter since the relative errors are sufficient.

Comparison of theoretical and numerical formulations
First a comparison between different mixed-formulation methods is presented, to give a reference and also to illustrate an improvement due to a suggested approach. In the subsequent figures, the eigenfrequencies are plotted within the complex plane where the real part {f } is associated with the corresponding physical eigenfrequency of the corresponding harmonic oscillation. Figures 7 and 8 show the results for three mixedformulations, namely solving the test case with LNSE, cf. Figure 7a, with LEE, cf. Figure 7b, and the mixed formulation of Galbrun's equation, cf. Figure 8. From Figure 7a and 7b, it is not possible to clearly identify the acoustic eigenfrequencies at f 1 = 50 Hz and f 2 = 100 Hz at Ma = 0 for both methods using LNSE and LEE. Moreover, the imaginary parts of the eigenvalues have small values leading to a possible pollution of the results when a subsequent modal superposition procedure is applied in order to restore the system response in the frame of a harmonic analysis. In contrast, the mixed formulation of Galbrun's equation shows a more stable behavior since the imaginary parts of all acoustically not relevant eigenvalues (i.e. convective and spurious modes) increase their value with increasing Mach number. As mentioned before, the authors propose the use of a discontinuous Galerkin discretization method to determine solutions to Galbrun's equation in a displacementbased formulation, rather than solving the mixed formulation with the use of continuous Galerkin discretization.
For this purpose, the two different discretization schemes are compared in Figure 9 and Figure 10 specifically to highlight the benefits of the DG method. In Figure 9a and 9b, the computed eigenvalues are plotted in the complex plane when using a continuous discretiza- tion scheme while utilizing linear and quadratic basis functions, respectively. It is observed that the eigenvalues are widely spread across the complex plane and the acoustic modes of interest at f 1 = 50 Hz and f 2 = 100 Hz are not distinguishable from the others.
In contrast, cf. Figure 10a and 10b, the use of a discontinuous Galerkin discretization method enables a clear separation of the relevant acoustic eigenfrequencies from the others.
Comparing the solution of the acoustic eigenfrequencies from the numerical examples of "GAL w1 DG" (linear elements) and "GAL w2 DG" (quadratic elements), cf. Figure 10a and 10b, with the analytical values according to equation (25), we find good agreement. These values are listed in Table II where the relative error is calculated in accordance to where f denotes the analytic result andf the numerical approximation.

Flux constant
As mentioned in the previous section, an appropriate value for α within the Lax-Friedichs-Flux needs to be chosen in order to adjust for a correct flux representation across each finite element. For this purpose, the influence of the parameter α with respect to the results is investigated for the filtered Galbrun equation, i.e. a Lagrange multiplier λ h is used to enforce ∇ × w = 0. By varying the parameter α for a constant Mach number, it is observed that for sufficiently accurate results, α should be larger than 10 6 , cf. Figure 11 and Figure 12, independent of flow velocity and order of basis function. The authors identified that α scales with (c 0 + v 0 ) 2 . These results are in close agreement to reported conclusions, see [13].

Vortical and spurious modes
The originality of Galbrun's equation is that the mixed pressure-displacement formulation given by Treyssède [41] is not changed by the presence of flow and is generally identical to the no-flow case. Except the presence of flow complicates the analysis to the existence of modal solutions due to the convective terms present. Even for an irrotational source term, the displacement field is not generally curl-free for shear-flow problems. Nevertheless, it has not yet been proven whether Discontinuous Galerkin methods satisfy the inf-sup condition. Vortical modes are non-acoustic perturbations convected with the mean flow, i.e. propagate with the flow. If the mean-flow is uniform, acoustic and vortical modes are decoupled but when the mean flow is rotational, the two types of modes couple. Hence it is possible to identify these modes from numerical artefacts, the spurious modes.
Further, especially for uniform mean-flow, it is clear that while the location of acoustic eigenfrequencies decrease with increasing Mach number, some other (non-acoustic) eigenfrequencies show a clear proportional dependency on the Mach number while others seem to be randomly distributed across the complex plane. Now, Figure 13 displays a similar version of Figure 10b. In Figure 13, the Mach number was sampled in 100 steps between Ma = 0 and Ma = 0.3, where each color code represents a different Mach number sample. This way, it is possible to follow the development of the eigenvalues with increasing Mach number. Figure 13  number, the eigenvector remains stable. The authors understand "stable" in the sense that the eigenvalue depends on the Mach number and is following a certain characteristic where the associated eigenvector remains unchanged in its appearance. The mentioned characteristics are frequency as well as mesh size dependent. Until today it remains an open question how these characteristics can be described and whether their nature is either physical or numerical? Notably, slightly altering the finite element mesh structure changes these Mach-number-dependent eigenvectors entirely.
To investigate this behavior further more, Figure 14 shows the eigenvectors of the two acoustic modes and two other modes that follow a certain characteristic Θ 3 (f 3 ) = const. and Θ 4 (f 4 ) = const., cf. Figure 14a. Since the acoustic pressure is proportional to the divergence of the displacement field, see equation (8), we present the divergence (∇ · w) and the rotation (∇ × w) of the displacement field in Figures 14b to 14i. It can be seen that for the acoustic modes, the rotation of the displacement field is two orders of magnitude lower compared to the divergence of the displacement field. In contrast, for non-acoustic modes the rotation is two orders of magnitude higher. This comparison is valid since the data to calculate the divergence and the rotation of the displacement field is relative to the mode.

Filtering Galbrun's equation
Bearing in mind results from the previous section, a filtering of Galbrun's equation is conducted to exclude the vorticity related eigenvalues from the solution space. Keeping this in mind, the Lagrange multiplier λ h formulation introduced, introduced in an earlier section, to enforce the rotation of the displacement field to be zero within the domain. The results are shown in Figure 15 and Figure 16.
It can be seen that by filtering Galbrun's equation, vorticity related eigenvalues are suppressed leaving only for acoustics relevant modes. Figures 15(b) and Figure 16(b) display an extended area of the complex plane. For larger imaginary values some eigenvalues persist. Since these eigenvalues are unlikely to be related to vorticity modes

Examples
In this section the capabilities of the numerical method to solve more realistic problems are illustrated. The first ex- ample represents a duct with uniform flow such as in Section 3 is considered. An admittance boundary condition is applied at the outlet of the duct. The remaining boundaries are considered as acoustically hard walls, i.e. the associated eigenvectors to the solution of Galbrun's equation have zero displacement on these boundaries. The second example can be seen as a cross section of a swirling flow represented by an annulus with a shear flow in circumferential direction. A boundary admittance is considered on the outer ring.

Finite duct with an absorbing end-condition
To extend our verification of the non-mixed Galbrun formulation, we present solutions for an equivalent boundary condition for Galbrun's equation in comparison to the Robin boundary condition when considering the Helmholtz equation. Figure 18 displays the model under consideration with the corresponding results. In this step the Mach number is set to zero in order to compare the results of the proposed method for solving Galbrun's equation with the standard Galerkin (CG) discretization of the Helmholtz equation using quadratic elements for the pressure unknowns. The mesh is chosen as illustrated in Figure 4.
As can be seen, the results are in very good agreement. Since the duct configuration is such that only plane wave propagation is possible below a cut on frequency of f = 340 Hz, there is no propagating mode below this frequency, which correlates with the results pictured above. For the acoustic eigenvalues inside a duct below the cuton frequency, Marburg [28] presented comparable results. In addition considering the dispersion relation, the same results, i.e. f = 340 Hz can be found.
When the flow velocity is increased, the eigenvalues spread into the complex plane, see Figure 19 number the convectional effect on the mode shape is increasing, too.

An annulus supporting a rotating shear flow
To investigate the capabilities of the proposed method on a more realistic problem, an annulus with circulating flow is analyzed. Figure 17 shows the configuration of the test case model. An annulus with an inner radius r 1 = 0.75 m and an outer radius of r 2 = 1 m is chosen, which is equivalent to the model investigated in Dietzsch et al. [12].
The mean flow velocity is defined as where r denotes the radial or "1"-direction. For simplicity the mean flow mass density is ρ 0 = 1 kg/m 3 and the mean pressure is set to p 0 = 0. For the case of Ma = 0.3 the mean flow would circulate clockwise at r = r 2 and counterclockwise at r = r 1 the tangential component of the mean flow velocity at r = r 2 would be |v 02 | = 102 m/s. The radial component remains zero at all the time. Again, the eigenvalues at four different flow velocities are investigated.
As a first case the boundary condition on the outer surface takesȲ = 0 which corresponds to an acoustically hard wall or w k n k = 0, respectively. In this case, the acoustically relevant modes can be identified with respect to their vanishing imaginary part. Figure 20 displays the results of the eigenvalues for the annulus configuration with acoustically hard walls. It is noticeable that all acoustic eigenvalues have a negligible imaginary part and can easily be separated. In addition, for Ma = 0, the eigenvalue at (f ) = 124.098 is a double mode due to the double symmetry of the geometry. As the mean flow velocity increases, the modes separate into a forward and a backward traveling wave where one eigenvalue increases and the other one decreases. The same behavior can be found when investigating tyre cavities under rotation, see Lopez et al. [26].
In a second example,Ȳ = 1 and the (∇ × w)filtered Galbrun equation is used. Figure 21 shows the results for eigenvalue computation. It can be seen that when increasingȲ , the eigenvalue distribution is shifted to lower imaginary values and is not symmetric with respect to the real axis. Further depending on the Mach number, characteristic lines are identifiable where eigenvalues group along. Investigating this further more, the authors identify acoustically relevant eigenvalues between these characteristic lines for a given Mach number. Additionally, for the acoustically relevant eigenvalues, a behavior where (∇ · w) (∇ × w) was noticed. For all other eigenvalues the divergence of the displacement field is in the same order of magnitude as the curl of the displacement field, i.e. (∇ · w) ≈ (∇ × w). This behavior can be used to separate acoustically relevant eigenvalues from others. At this point it must be noted that with a pure displacement based formulation of Galbrun's equation such a determination is easily possible. So any scaling of the eigenvectors would result in a scaling of the divergence and the curl of the displacement field.

Conclusions
When solving Galbrun's equation in its pure displacement based formulation, the authors propose the use of a discontinuous Galerkin method for discretizing Galbrun's equation. Due to discontinuity between neighboring elements it is shown that a Lax-Friedrichs condition with a flux factor α ≥ 10 6 gives sufficiently accurate results. In order to exclude vorticity modes from the modal space, a (∇ × w)filtering with the aid of a Lagrange multiplier is successfully applied and illustrating by various examples. Since the pressure distribution is related to the divergence of the displacement field, the presented studies show the separation between acoustic modes from vorticity modes or numerical spurious modes. It is shown that the vorticity modes follow a certain characteristic which are frequency  and mesh size dependent. To account for admittance conditions on domain boundaries, an adequate boundary formulation is presented.
Future research will be dedicated towards formulating appropriate boundary conditions for solving general exterior problems and duct problems by robust reformulations provided by Bonnet-Benn Dhia et al. [5], for example, using discontinuous Galerkin numerical solutions.