2.1. Direct problem

2.1.1. Governing equation

Goldstein (Reference Goldstein1978) derived the following scalar wave equation for the perturbation $\phi$ of the velocity potential:

(2.1)\begin{equation} \rho_0 {\frac{\mathrm{d}_0}{\mathrm{d} t}} \left( \frac{1}{c_0^2}\,{\frac{\mathrm{d}_0 \phi}{\mathrm{d} t}} \right) - \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \rho_0\,\boldsymbol{\nabla}\phi \right) ={-}q , \end{equation}

where ${\mathrm {d}_0/\mathrm {d} t}=\partial /\partial t + \boldsymbol {u}_0 \boldsymbol {\cdot } \boldsymbol {\nabla }$ is the material derivative with respect to the mean flow velocity $\boldsymbol {u}_0$, $c_0$ is the speed of sound, $\rho _0$ is the mean flow density, and $q$ is a generic distributed source. From the potential $\phi$, it is possible to compute the other acoustic quantities such as pressure, density and velocity:

(2.2a–c)\begin{equation} p ={-}\rho_0\,{\frac{\mathrm{d}_0\phi}{\mathrm{d} t}} ,\quad \rho=p/c_0^2,\quad \boldsymbol{u}=\boldsymbol{\nabla}\phi. \end{equation}

This model is solved in the frequency domain using the implicit time dependence $\mathrm {e}^{+\mathrm {i}\omega t}$. This amounts to replacing the material derivative $\mathrm {d}_0/\mathrm {d} t$ in the above expressions by its frequency-domain counterpart $\mathrm {D}_0/\mathrm {D} t=\mathrm {i}\omega + \boldsymbol {u}_0 \boldsymbol {\cdot } \boldsymbol {\nabla }$.

2.1.2. Variational formulation

For a domain $\varOmega$ with boundary $\partial \varOmega$, the variational formulation for (2.1) reads

(2.3)\begin{align} &\int_{\varOmega} \left(\rho_0\,\boldsymbol{\nabla} \bar{\psi} \boldsymbol{\cdot} \boldsymbol{\nabla} \phi -\frac{\rho_0}{c_0^2}\,\overline{\frac{\mathrm{D}_0 \psi}{\mathrm{D} t}}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t} \right) \mathrm{d} \varOmega + \int_{\partial\varOmega} \rho_0\bar{\psi}\left(\frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{n}} {c_0^2}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t} - \frac{\partial \phi}{\partial n} \right)\mathrm{d} \varGamma \nonumber\\ &\quad ={-}\int_{\varOmega} \bar{\psi} q \,\mathrm{d} \varOmega , \end{align}

where $\psi$ is the test function associated with the velocity potential $\phi$, $\bar {\cdot }$ is the complex conjugate, and $\boldsymbol {n}$ is the unit outward normal vector on $\partial \varOmega$.

The boundary integral in (2.3) should be investigated to consider each different boundary condition. In this paper, we will consider two types of boundary conditions: a lined and vibrating surface denoted $\varGamma$, and a radiation condition imposed on a boundary $\varGamma _\infty$ located in the far field. We now describe in detail these two boundary conditions.

2.1.3. Boundary condition for a lined and vibrating surface

The boundary $\varGamma$ is a lined and vibrating surface with acoustic admittance $A$ and a prescribed normal velocity $v_s$. It is also assumed impervious to the mean flow, therefore $\boldsymbol {u}_0\boldsymbol {\cdot }\boldsymbol {n}=0$. The corresponding boundary integral in the variational formulation (2.3) is

(2.4)\begin{equation} I_{\varGamma} ={-}\int_{\varGamma} \rho_0 \bar{\psi}\,\frac{\partial \phi}{\partial n}\,\mathrm{d} \varGamma . \end{equation}

Due to the grazing mean flow, the presence of the acoustic treatment on the surface is accounted for with the Myers (Reference Myers1980) condition, which assumes an infinitely thin boundary layer above the surface. It relates the normal fluctuating velocity $\boldsymbol {u}\boldsymbol {\cdot }\boldsymbol {n}$ of the fluid above the boundary layer to the normal displacement $\xi$ of the fluid on the surface:

(2.5)\begin{equation} \boldsymbol{u}\boldsymbol{\cdot}\boldsymbol{n} = \frac{\mathrm{D}_0\xi}{\mathrm{D} t} - \xi \boldsymbol{n}\boldsymbol{\cdot} [(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla} )\boldsymbol{u}_0], \end{equation}

where $\mathrm {i}\omega \xi$ is the total normal velocity of the surface defined as a combination of the prescribed normal velocity $v_s$ and the relative fluid velocity ${A} p$ through the surface allowed by the acoustic treatment, i.e.

(2.6)\begin{equation} \mathrm{i}\omega\xi = v_s +{A} p . \end{equation}

Upon introducing the Myers condition (2.5), the boundary term (2.4) becomes

(2.7)\begin{equation} I_{\varGamma} ={-} \int_{\varGamma} \rho_0 \bar{\psi} \left\{ \frac{\mathrm{D}_0\xi}{\mathrm{D} t} - \xi \boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n} \boldsymbol{\cdot} \boldsymbol{\nabla})\boldsymbol{u}_0] \right\}\mathrm{d} \varGamma. \end{equation}

Eversman (Reference Eversman2001a) shows that this integral can be simplified by using the following result from vector analysis (see also Möhring Reference Möhring2001),

(2.8)\begin{equation} \rho_0 \bar{\psi} \xi \boldsymbol{n}\boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0] = \boldsymbol{\nabla}\boldsymbol{\cdot}(\rho_0\bar{\psi}\xi\boldsymbol{u}_0) - \boldsymbol{n}\boldsymbol{\cdot}[\boldsymbol{\nabla}\times (\boldsymbol{n}\times\rho_0\bar{\psi} \xi \boldsymbol{u}_0)] , \end{equation}

together with Stokes’ theorem. This yields the following expression for the boundary integral:

(2.9)\begin{equation} I_\varGamma = \int_{\varGamma} \rho_0\,\overline{\frac{\mathrm{D}_0\psi}{\mathrm{D} t}}\,\xi \,\mathrm{d} \varGamma + \oint_{\partial\varGamma} \rho_0\bar{\psi}\xi(\boldsymbol{u_0}\times\boldsymbol{n})\boldsymbol{\cdot}\mathrm{d}\boldsymbol{l} ,\end{equation}

which involves a line integral along the contour of the boundary $\varGamma$. Following Eversman (Reference Eversman2001a) and Rienstra (Reference Rienstra2007), it can be argued that the normal displacement $\xi$ should be continuous between a lined surface and a rigid surface (to ensure the continuity of the unsteady streamlines). The displacement should therefore vanish on the contour $\partial \varGamma$, and the contour integral can be removed. In the present case, this implies that both the admittance $A$ and the velocity $v_s$ should vanish on the contour of the surface $\varGamma$. With these assumptions, the contour integral along $\partial \varGamma$ is dropped in what follows. Note that this has an impact on the solution, mostly in the vicinity of the liner discontinuity and on the reflection coefficient (Gabard Reference Gabard2010).

After using (2.5) and (2.2a–c), the boundary integral for a lined vibrating surface finally reads

(2.10)\begin{equation} I_{\varGamma} = \int_{\varGamma} \frac{\rho_0}{\mathrm{i}\omega}\, \overline{\frac{\mathrm{D}_0 \psi}{\mathrm{D} t}}\,v_s \,\mathrm{d} \varGamma -\int_{\varGamma} \frac{\rho_0^2 {A}}{\mathrm{i} \omega }\,\overline{\frac{\mathrm{D}_0\psi}{\mathrm{D} t}}\,\frac{\mathrm{D}_0\phi}{\mathrm{D} t} \,\mathrm{d} \varGamma . \end{equation}

The first integral is a forcing term that will appear on the right-hand side of the discretised system.

2.1.4. Far-field radiation condition

We now consider an outer surface $\varGamma _{\infty }$ located far away from the other surfaces and the volume sources. The contribution from this surface to the variational formulation is

(2.11)\begin{equation} I_{\varGamma_{\infty}} = \int_{\varGamma_{\infty}} \rho_0\bar{\psi}\left(\frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{n}} {c_0^2}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t} - \frac{\partial \phi}{\partial n} \right)\mathrm{d} \varGamma . \end{equation}

Without loss of generality, $\varGamma _\infty$ is chosen as the sphere defined by $\|\boldsymbol {x}\|=r$, hence $\boldsymbol {n}=\boldsymbol {x}/r$ and $\partial /\partial n=\partial /\partial r$.

For the Helmholtz equation (i.e. with no mean flow), one is left with $\partial \phi /\partial r$ on this outer boundary, and the Sommerfeld radiation condition states that the asymptotic behaviour of $\phi$ in the far field is such that $\partial \phi /\partial r = -\mathrm {i} k \phi$ with $k=\omega /c_0$. The Sommerfeld radiation condition can be generalised to include the effect of a mean flow (Bayliss & Turkel Reference Bayliss and Turkel1982; Tam & Webb Reference Tam and Webb1993; Bogey & Bailly Reference Bogey and Bailly2002). However these radiation conditions are not directly applicable to the variational formulation (2.3). To derive the radiation condition in a form suitable for (2.1), we have to determine the asymptotic behaviour of all the terms in the parentheses in (2.11).

To that end, we consider the sound field radiated in a uniform mean flow by a generic source term $q(\boldsymbol {x})$ on the right-hand side of (2.1). Using Green's formula, the radiated sound field can be written

(2.12)\begin{equation} \phi(\boldsymbol{x}) ={-}\int q(\boldsymbol{y})\,G(\boldsymbol{x}|\boldsymbol{y})\,\mathrm{d}\boldsymbol{y} , \end{equation}

where the free-field Green's function $G(\boldsymbol {x}|\boldsymbol {y})$ is defined in three dimensions by

(2.13)\begin{equation} G(\boldsymbol{x}|\boldsymbol{y}) = \frac{1}{4 {\rm \pi}\tilde{r}}\exp\left\{-\frac{\mathrm{i} k}{\beta^2} \left[\tilde{r} - (\boldsymbol{x}-\boldsymbol{y})\boldsymbol{\cdot}\frac{\boldsymbol{u}_0}{c_0}\right] \right\} , \end{equation}

with $\beta ^2=1-\|\boldsymbol {u}_0\|^2/c_0^2$, and the distance $\tilde {r}$ defined as

(2.14)\begin{equation} \tilde{r} = \sqrt{\beta^2 \|\boldsymbol{x} - \boldsymbol{y}\|^2 + [(\boldsymbol{x}-\boldsymbol{y})\boldsymbol{\cdot}\boldsymbol{u}_0]^2/c_0^2} . \end{equation}

Note that this Green's function already satisfies the far-field radiation condition. After some lengthy developments, it is possible to obtain the following result:

(2.15)\begin{equation} \frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{n}}{c_0^2}\,\frac{\mathrm{D}_0\,G(\boldsymbol{x}|\boldsymbol{y})}{\mathrm{D} t} - \frac{\partial G(\boldsymbol{x}|\boldsymbol{y})}{\partial n} = \mathrm{i} k\,\frac{r}{\tilde{r}} \left( 1 + \frac{\beta^2}{\mathrm{i} k\tilde{r}} \right) \left( 1 - \frac{\boldsymbol{x}\boldsymbol{\cdot}\boldsymbol{y}}{r\tilde{r}} \right) G(\boldsymbol{x}|\boldsymbol{y}) , \end{equation}

where the derivatives on the left-hand side operate on the $\boldsymbol {x}$ coordinate. Note that this result is exact and valid for any value of $r$, i.e. no far-field approximation has been made.

When the observer is in the geometric far field, which is defined by $\|\boldsymbol {x}\|\gg \|\boldsymbol {y}\|$, the term $\boldsymbol {x}\boldsymbol {\cdot }\boldsymbol {y}/(r\tilde {r})$ in (2.15) can be neglected, and $\tilde {r} \simeq \hat {r} = \sqrt {\beta ^2 \|\boldsymbol {x}\|^2 + (\boldsymbol {x}\boldsymbol {\cdot }\boldsymbol {u}_0)^2/c_0^2}$. From (2.12) and (2.15), it is apparent that the velocity potential satisfies the following radiation condition in the geometric far field:

(2.16)\begin{equation} \frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{x}}{c_0^2 r}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t} - \frac{\partial \phi}{\partial r} = \mathrm{i} k\,\frac{r}{\hat{r}} \left( 1 + \frac{\beta^2}{\mathrm{i} k\hat{r}} \right) \phi + {O}\left( \frac{1}{r^2}\right) ,\quad\text{when }r\to+\infty . \end{equation}

In the absence of mean flow, this radiation condition reduces to the standard Sommerfeld condition; see, for instance, § 4.5.4 in Pierce (Reference Pierce2019). The corresponding boundary integral in the variational formulation is written

(2.17)\begin{equation} I_{\varGamma_\infty} = \int_{\varGamma_\infty}\mathrm{i} k\alpha\rho_0\bar{\psi}\phi\,\mathrm{d}\varGamma , \end{equation}

where the coefficient $\alpha$ is defined as

(2.18)\begin{equation} \alpha = \frac{r}{\hat{r}} \left( 1 + \frac{\beta^2}{\mathrm{i} k\hat{r}} \right) . \end{equation}

In (2.16), we have kept the leading-order term, but it is also possible to derive higher-order radiation conditions from (2.15); see Givoli (Reference Givoli2004).

The corresponding radiation condition for a two-dimensional problem is given in Appendix A.

2.2. The adjoint operator and tailored Green's functions

As explained above, it can be particularly useful to use a tailored adjoint Green's function to facilitate or accelerate the calculation of solutions to the direct problems, in this case by providing an explicit expression for the acoustic potential $\phi$ in terms of the sources $v_s$ and $q$. Such a tailored Green's function is a solution to an adjoint problem that will be devised in this subsection.

We first introduce scalar products between complex-valued functions, on either $\varOmega$ or $\varGamma$, as follows:

(2.19a,b)\begin{equation} \left\langle f,g\right\rangle_{\varOmega}=\int_{\varOmega} \overline{f(\boldsymbol{y})}\,g(\boldsymbol{y}) \,\mathrm{d} \varOmega_y \quad \left\langle f,g\right\rangle_{\varGamma}= \int_{\varGamma} \overline{f(\boldsymbol{y})}\, g(\boldsymbol{y}) \,\mathrm{d} \varGamma_y . \end{equation}

We then form the scalar product of Goldstein's equation (2.1) with a generic Green's function $G_{\phi }(\boldsymbol {x},\boldsymbol {y})$, which remains to be defined:

(2.20)\begin{align} &\int_{\varOmega} \overline{G_{\phi}}(\boldsymbol{x}, \boldsymbol{y})\, \rho_0(\boldsymbol{y})\,\frac{\mathrm{D}_0}{\mathrm{D} t}\left[\frac{1}{ c_0^2(\boldsymbol{y}) }\, \frac{\mathrm{D}_0}{\mathrm{D} t}\phi(\boldsymbol{y}) \right] -\overline{G_{\phi}}(\boldsymbol{x}, \boldsymbol{y})\,\boldsymbol{\nabla} \boldsymbol{\cdot} \left[\rho_0(\boldsymbol{y})\,\boldsymbol{\nabla} \phi(\boldsymbol{y}) \right]\mathrm{d}{\varOmega_y}\nonumber\\ &\quad ={-}\int_{\varOmega}\overline{G_{\phi}}(\boldsymbol{x}, \boldsymbol{y})\,q(\boldsymbol{y}) \,\mathrm{d}{\varOmega_y} , \end{align}

which is valid for any point $\boldsymbol {x}$ located in $\varOmega$, where the integrals are performed over $\boldsymbol {y}$. In the following expressions, we will omit the $\boldsymbol {x}$ and $\boldsymbol {y}$ dependence to simplify the notation. However, to avoid ambiguities, note that all quantities are functions of only $\boldsymbol {y}$ except for $G_\phi$, and all the derivatives are applied with respect to $\boldsymbol {y}$. Integrating twice by parts the expression above yields

(2.21)\begin{align} &\int_{\varOmega} \phi \overline{\left[\rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{ c_0^2}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} \right) -\boldsymbol{\nabla}_y \boldsymbol{\cdot} \left( \rho_0\,\boldsymbol{\nabla}_y G_{\phi}\right) \right]}\,\mathrm{d} \varOmega_y\nonumber\\ &\qquad +\int_{\partial\varOmega} \frac{\rho_0(\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{n})}{c_0^2} \left( \overline{G_{\phi}}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t}- \phi\, \overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}} \right) + \rho_0 \left( \phi\, \frac{\partial \overline{G_{\phi}}}{\partial {n_y}} -\overline{G_{\phi}}\, \frac{\partial\phi}{\partial {n_y}} \right)\mathrm{d} \varGamma_y\nonumber\\ &\quad ={-}\int_{\varOmega}\overline{G_{\phi}} q \,\mathrm{d}\varOmega_y . \end{align}

Subtracting (2.21) and (2.20) yields the integral form of Lagrange's identity, which relates the direct and adjoint operators; see, for instance, § IV.4 in Dennery & Krzywicki (Reference Dennery and Krzywicki2012). The adjoint operator to Goldstein's equation (2.1) is readily found in the square brackets in (2.21). In this case, the direct and adjoint operators are identical, which is expected since Goldstein's equation is self-adjoint.

2.2.1. Adjoint equation

Our aim is now to identify a tailored adjoint Green's function that provides an explicit solution for the direct problem defined in the previous subsection. In other words, we have to select the governing equation and the boundary conditions for $G_\phi$ so that the velocity potential $\phi$ can be written explicitly in terms of the source terms $q$ and $v_s$. Beginning with the governing equation for $G_\phi$, we would like the first integral in (2.21) to reduce to $\phi (\boldsymbol {x})$, i.e. $G_\phi$ should be such that

(2.22)\begin{equation} \int_{\varOmega} \phi \overline{\left[\rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{ c_0^2}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} \right) -\boldsymbol{\nabla}_y \boldsymbol{\cdot} \left( \rho_0\,\boldsymbol{\nabla}_y G_{\phi}\right) \right]}\,\mathrm{d} \varOmega_y = \phi(\boldsymbol{x}) . \end{equation}

This can be achieved if $G_\phi$ is a solution of

(2.23)\begin{equation} \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{ c_0^2}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} \right) -\boldsymbol{\nabla}_y \boldsymbol{\cdot} \left( \rho_0\,\boldsymbol{\nabla}_y G_{\phi}\right) = \delta(\boldsymbol{x}-\boldsymbol{y}), \end{equation}

which shows that $G_\phi$ is indeed a Green's function for the adjoint operator.

Looking back at (2.21), the solution $\phi (\boldsymbol {x})$ is now given in terms of integrals on $\varGamma$ and $\varOmega$. While the latter depends on only $q$ and $G_\phi$, the boundary integral still involves $\phi$. In order to obtain an explicit expression for $\phi$, we have to choose the boundary conditions for the Green's function so that the boundary integral involves only $G_\phi$ and $v_s$.

2.2.2. Lined vibrating wall

We first consider the integral on the boundary $\varGamma$ representing a lined and vibrating surface, as described in § 2.1.3. Since $\boldsymbol {u}_0\boldsymbol {\cdot }\boldsymbol {n}=0$ on this surface, the relevant boundary integral in (2.21) is

(2.24)\begin{equation} J_{\varGamma} =\int_{\varGamma } \rho_0 \left( \phi\,\frac{\partial \overline{G_{\phi}}}{\partial {n_y}} -\overline{G_{\phi}}\, \frac{\partial\phi}{\partial {n_y}} \right)\mathrm{d} \varGamma_y. \end{equation}

The second term in this integral is the same as in (2.4) with $\psi$ replaced by $G_\phi$. Using the same analysis as in § 2.1.3 for the formulation of the Myers boundary condition, one gets

(2.25)\begin{equation} J_{\varGamma} = \int_\varGamma \frac{\rho_0}{\mathrm{i} \omega}\, \overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,v_s + \rho_0\phi\,\frac{\partial \overline{G_{\phi}}}{\partial {n_y}} - \frac{\rho_0^2{A}}{\mathrm{i} \omega }\,\overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,\frac{\mathrm{D}_0 \phi}{\mathrm{D} t} \,\mathrm{d} \varGamma . \end{equation}

To proceed further, we use steps similar to those described in § 2.1.3 for the formulation of the Myers condition: i.e. the derivative in $\mathrm {D}_0\phi /\mathrm {D} t$ is integrated by parts, an expression similar to (2.8) is used, and we apply Stokes’ theorem. This yields

(2.26)\begin{align} J_{\varGamma} &= \int_\varGamma \frac{\rho_0}{\mathrm{i} \omega}\, \overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,v_s \,\mathrm{d} \varGamma +\oint_{\partial\varGamma} \phi\,\frac{\rho_0^2{A}}{\mathrm{i} \omega}\,\overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,(\boldsymbol{u}_0\times\boldsymbol{n})\boldsymbol{\cdot}\mathrm{d}\boldsymbol{l}\nonumber\\ &\quad +\int_{\varGamma} \rho_0 \phi\,\overline{\left\{ \frac{\partial G_{\phi}}{\partial {n_y}} - \frac{\mathrm{D}_0 }{\mathrm{D} t}\left(\frac{\rho_0\bar{{A}}}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}\right) + \frac{\rho_0\bar{{A}}}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}\,\boldsymbol{n}\boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot} \boldsymbol{\nabla}_y)\boldsymbol{u}_0] \right\}} . \end{align}

Recall that the aim is to remove $\phi$ from this expression. First, we can use the same argument as in Eversman (Reference Eversman2001a) for the direct problem: assuming that the admittance $A$ varies continuously on $\varGamma$, the contour $\partial \varGamma$ can be located on the rigid wall and the contour integral above vanishes. Second, the last integral in (2.26) is eliminated by choosing the following boundary condition for the tailored adjoint Green's function $G_\phi$:

(2.27)\begin{equation} \frac{\partial {G_{\phi}}}{\partial {n_y}} = \frac{\mathrm{D}_0\eta}{\mathrm{D} t} - \eta\boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot} \boldsymbol{\nabla}_y) \boldsymbol{u}_0], \quad \text{with } \mathrm{i}\omega\eta =\rho_0 \bar{{A}}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} ,\quad\text{on }\varGamma. \end{equation}

This expression is the adjoint Myers condition. It takes a form similar to the original impedance condition (2.5) for the direct problem except for several important differences. There is a sign difference compared to (2.5), and the admittance $A$ is replaced by its complex conjugate. In addition, the source term $v_s$ does not appear in the adjoint impedance condition (2.27) for the tailored Green's function $G_\phi$.

2.2.3. Radiation condition

The contribution of the radiation boundary $\varGamma _\infty$ in (2.21) is

(2.28)\begin{equation} J_{\varGamma_\infty} = \int_{ \varGamma_{\infty}} \frac{\rho_0(\boldsymbol{u}_0\boldsymbol{\cdot} \boldsymbol{n})}{c_0^2} \left( \overline{G_{\phi}}\, \frac{\mathrm{D}_0 \phi}{\mathrm{D} t}- \phi\,\overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}} \right) + \rho_0 \left( \phi\,\frac{\partial \overline{G_{\phi}}}{\partial {n_y}} -\overline{G_{\phi}}\,\frac{\partial\phi}{\partial {n_y}} \right)\mathrm{d} \varGamma_y . \end{equation}

Using the radiation condition (2.16) satisfied by $\phi$ on $\varGamma _{\infty }$, this contribution becomes

(2.29)\begin{equation} J_{\varGamma_\infty} = \int_{ \varGamma_{\infty}} \rho_0 \phi \overline{\left( \frac{\partial G_{\phi}}{\partial {n_y}} - \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{n}}{c_0^2}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} -\mathrm{i} k\overline{\alpha} G_{\phi} \right)}\, \mathrm{d} \varGamma_y , \end{equation}

with $\alpha$ defined in (2.17). For this contribution to vanish, the Green's function $G_\phi$ should satisfy the following adjoint radiation condition on $\varGamma _{\infty }$:

(2.30)\begin{equation} \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{y}}{c_0^2 r}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} - \frac{\partial G_{\phi}}{\partial r} ={-}\mathrm{i} k \bar{\alpha} G_{\phi} . \end{equation}

It is similar to the direct radiation condition (2.16) except that $\mathrm {i} k \alpha$ is replaced by its complex conjugate $-\mathrm {i} k \bar {\alpha }$. This is a consequence of the fact that the adjoint problem is anti-causal. The adjoint radiation condition allows only inward-propagating waves, which is the opposite of the radiation condition for the direct problem, which allows only outward-propagating waves.

2.2.4. Summary and discussion

The acoustic potential $\phi$ can be written explicitly in terms of the source terms $v_s$ and $q$ as

(2.31)\begin{equation} \phi(\boldsymbol{x}) = \left\langle -\rho_0\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t},\frac{v_s}{\mathrm{i} \omega} \right\rangle_{\varGamma} - \left\langle G_{\phi},q\right\rangle_{\varOmega} ={-}\int_{\varGamma} \frac{\rho_0}{\mathrm{i} \omega}\,\overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,v_s\,\mathrm{d}{\varGamma_y} -\int_{\varOmega}\overline{G_{\phi}} q \,\mathrm{d}{\varOmega_y} , \end{equation}

provided that the tailored adjoint Green's function $G_{\phi }$ is defined as follows:

(2.32)\begin{equation} \left.\begin{aligned} & \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t}\left(\frac{1}{ c_0^2 }\,\frac{\mathrm{D}_0}{\mathrm{D} t}G_{\phi}\right) - \boldsymbol{\nabla}_y \boldsymbol{\cdot} \left(\rho_0\,\boldsymbol{\nabla}_y G_{\phi} \right) = \delta(\boldsymbol{x} -\boldsymbol{y}) \quad\text{in }\varOmega,\\ & \frac{\partial {G_{\phi}}}{\partial {n_y}} = \frac{\mathrm{D}_0\eta}{\mathrm{D} t} - \eta\boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot} \boldsymbol{\nabla}_y) \boldsymbol{u}_0] \quad \text{with } \mathrm{i}\omega\eta =\rho_0 \bar{{A}}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} \quad\text{on } \varGamma,\\ & \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{y}}{c_0^2 r}\,\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t} - \frac{\partial G_{\phi}}{\partial r} ={-}\mathrm{i} k \bar{\alpha} G_{\phi} \quad \text{on }\varGamma_\infty. \end{aligned} \right\} \end{equation}

It is worth noting a number of points concerning the direct and adjoint problems.

1. (i) While the propagation operator is self-adjoint, the boundary conditions are not self-adjoint since they differ between the direct and adjoint problems. The practical implications of these differences will be discussed in § 4.

2. (ii) The surface admittance $A$ does not appear explicitly in (2.31). It is taken into account implicitly by the tailored Green's function $G_{\phi }$ through the adjoint impedance condition (2.27) that it satisfies on the surface $\varGamma$.

3. (iii) In the expression (2.31), it is clear that $\boldsymbol {x}$ represents the observer position while the source position $\boldsymbol {y}$ moves on the surface $\varGamma$ and in the volume $\varOmega$. However, in the definition (2.32) of the tailored adjoint Green's function, the differential equation is written in terms of $\boldsymbol {y}$ (the derivative in the equation and the boundary conditions are applied with respect to $\boldsymbol {y}$ and not $\boldsymbol {x}$). This means that in (2.32), the source position is $\boldsymbol {x}$ and the observer position is $\boldsymbol {y}$, which is the reverse from (2.31). The reciprocity principle appears naturally from the derivation of $G_{\phi }$ presented above.

4. (iv) The tailored adjoint Green's function $G_{\phi }$ in (2.31) can be understood as the acoustic transfer function between the volume source $q$ and the velocity potential $\phi$ at the observer position $\boldsymbol {x}$. Likewise, the quantity $-\rho _0\,\mathrm {D}_0G_\phi /\mathrm {D} t$ is the transfer function between the surface displacement $v_s/(\mathrm {i}\omega )$ and the velocity potential $\phi (\boldsymbol {x})$.

In the analysis above, the tailored adjoint Green's function $G_\phi$ was devised to yield an explicit expression for the velocity potential. If one wishes to compute the acoustic pressure instead, then it is possible to compute $p(\boldsymbol {x})=-\rho _0\,\mathrm {D}_0\phi /\mathrm {D} t$ (with derivatives with respect to $\boldsymbol {x}$) from (2.31). This can be cumbersome as (2.32) yields $G_\phi$ as a function of $\boldsymbol {y}$ for a fixed $\boldsymbol {x}$. Instead, it is preferable to define another tailored adjoint Green's function $G_p(\boldsymbol {x},\boldsymbol {y})$ that relates directly $p$ to the source terms $q$ and $v_s$ through

(2.33)\begin{equation} p(\boldsymbol{x}) = \left\langle-\rho_0\,\frac{\mathrm{D}_0 G_{p}}{\mathrm{D} t},\frac{v_s}{\mathrm{i} \omega} \right\rangle_{\varGamma} - \left\langle G_{p},q\right\rangle_{\varOmega} , \end{equation}

which has the exact same form as (2.31). This is achieved by modifying the right-hand side of (2.23) so that the integral in (2.22) yields $p(\boldsymbol {x})$; that is, the adjoint equation for $G_p$ reads

(2.34)\begin{equation} \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{ c_0^2}\,\frac{\mathrm{D}_0 G_{p}}{\mathrm{D} t} \right) -\boldsymbol{\nabla}_y \boldsymbol{\cdot} \left( \rho_0\,\boldsymbol{\nabla}_y G_{p}\right) = \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t}\delta(\boldsymbol{x}-\boldsymbol{y}) . \end{equation}

It is straightforward to show that $G_p$ satisfies the same adjoint boundary conditions as in (2.32). A similar methodology can be used to compute the acoustic velocity field $\boldsymbol {u}$ instead of $\phi (\boldsymbol {x})$.

A final comment is that the generic source term $q$ can be replaced easily by a distribution of dipole sources by writing $q=\boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {f}$, where $\boldsymbol {f}$ is the dipole strength. The additional divergence in (2.31) can be transferred onto the Green's function by integrating by parts to write

(2.35)\begin{equation} \phi(\boldsymbol{x}) ={-}\int_{\varGamma} \frac{\rho_0}{\mathrm{i} \omega}\, \overline{\frac{\mathrm{D}_0 G_{\phi}}{\mathrm{D} t}}\,v_s\,\mathrm{d}{\varGamma_y} +\int_{\varOmega}\boldsymbol{f}\boldsymbol{\cdot}\boldsymbol{\nabla}_y\,\overline{G_{\phi}} \,\mathrm{d}{\varOmega_y} . \end{equation}

A similar approach can be used for other types of sources, such as a quadrupole distribution.

3.1. Direct problem

Instead of the velocity potential $\phi$, Pierce (Reference Pierce1990) proposes to use the momentum potential $\varphi$, which is such that $p=\mathrm {D}_0\varphi /\mathrm {D} t$ and $\rho _0\boldsymbol {u}=-\boldsymbol {\nabla }\varphi$. In the high-frequency limit, this variable satisfies the propagation equation

(3.1)\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} \left(\frac{1}{\rho_0}\,\boldsymbol{\nabla} \varphi \right) -\rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t}\left(\frac{1}{\rho_0^2 c_0^2 }\, \frac{\mathrm{D}_0\varphi}{\mathrm{D} t} \right)={-} q , \end{equation}

which is valid for an arbitrary steady base flow, i.e. compressible, rotational and not necessarily homentropic. Compared to (27) in Pierce (Reference Pierce1990) we have added a generic source distribution on the right-hand side. For a computational domain $\varOmega$ with boundary $\partial \varOmega$, the variational formulation corresponding to (3.1) reads

(3.2)\begin{align} &\int_{\varOmega} \frac{1}{\rho_0}\left(\frac{1}{c_0^2}\,\overline{\frac{ \mathrm{D}_0 \psi}{\mathrm{D} t}}\,\frac{\mathrm{D}_0 \varphi}{\mathrm{D} t} - \boldsymbol{\nabla} \bar{\psi}\boldsymbol{\cdot}\boldsymbol{\nabla} \varphi \right) \mathrm{d} \varOmega + \int_{\partial\varOmega} \frac{\bar{\psi}}{\rho_0} \left( \frac{\partial\varphi}{\partial n} - \frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{n}}{c_0^2}\,\frac{D_0\varphi}{\mathrm{D} t} \right) \mathrm{d} \varGamma \nonumber\\ &\quad ={-}\int_{\varOmega}\bar{\psi} q \,\mathrm{d}\varOmega , \end{align}

with the test function $\psi$. Like in (2.3) for Goldstein's equation, the volume integral above is Hermitian, which is consistent with the fact that there is an energy conservation principle associated with Pierce's equation (3.1); see § 2.4 in Möhring (Reference Möhring1999).

The Myers condition is traditionally used with Goldstein's equation (2.1), but we show here how it can be included in the variational formulation for Pierce's equation. When formulated in terms of the potential $\varphi$, the Myers impedance condition (2.5) becomes

(3.3)\begin{equation} \frac{\partial\varphi}{\partial n} ={-}\rho_0 \left\{\frac{\mathrm{D}_0\xi}{\mathrm{D} t} - \xi\boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0]\right\} , \end{equation}

with

(3.4)\begin{equation} \mathrm{i}\omega\xi = v_s + {A}\,\frac{\mathrm{D}_0 \varphi}{\mathrm{D} t}, \end{equation}

where $\xi$ is again the normal displacement of the surface. Recalling that $\boldsymbol {u}_0\boldsymbol {\cdot }\boldsymbol {n}=0$ on $\varGamma$, the boundary integral on surface $\varGamma$ in (3.2) reduces to

(3.5)\begin{equation} I_{\varGamma} = \int_{\varGamma}\frac{\bar{\psi}}{\rho_0}\,\frac{\partial\varphi}{\partial n}\,\mathrm{d} \varGamma ={-}\int_{\varGamma}\rho_0 \left(\frac{\bar{\psi}}{\rho_0}\right)\left\{\frac{\mathrm{D}_0\xi}{\mathrm{D} t} - \xi\boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0]\right\}\mathrm{d} \varGamma . \end{equation}

Following the same analysis as in § 2.1.3 (in particular using (2.8) where $\psi$ is replaced by $\psi /\rho _0$), it is possible to write

(3.6)\begin{equation} I_{\varGamma} = \int_{\varGamma}\frac{\rho_0}{\mathrm{i} \omega}\,\overline{\frac{\mathrm{D}_0}{\mathrm{D} t}\left(\frac{\psi}{\rho_0}\right)}\,v_s \,\mathrm{d} \varGamma +\int_{\varGamma} \frac{\rho_0 {A}}{\mathrm{i} \omega }\, \overline{\frac{\mathrm{D}}{\mathrm{D} t}\left( \frac{\psi}{\rho_0}\right)}\,\frac{\mathrm{D} \varphi}{\mathrm{D} t} \,\mathrm{d} \varGamma . \end{equation}

This formulation of the Myers impedance condition for Pierce's propagation equation is consistent with Eversman (Reference Eversman2001a), i.e. the displacement $\xi$ and admittance ${A}$ are assumed to vary smoothly between a lined and rigid surfaces. The main difference with (2.10) is the substitution of $\psi$ by $\psi /\rho _0$, which will be discussed further in § 4.

Concerning the radiation boundary $\varGamma _{\infty }$, the analysis of the far-field behaviour of the velocity potential $\phi$ given in § 2.1.4 also applies to the potential $\varphi$. It therefore satisfies the same radiation condition (2.16) as $\phi$, and the boundary integral on $\varGamma _{\infty }$ in (3.2) becomes

(3.7)\begin{equation} I_{\varGamma_{\infty}} ={-}\int_{\varGamma_{\infty}} \mathrm{i} k\alpha\,\frac{\bar{\psi}\varphi}{\rho_0} \,\mathrm{d} \varGamma . \end{equation}

3.2. Adjoint problem

Pierce's equation (3.1) can also be solved using a tailored adjoint Green's function, following the same approach as described in § 2.2 for Goldstein's equation. One starts by forming the scalar product between Pierce's equation (3.1) and a Green's function $H_{\varphi }(\boldsymbol {x},\boldsymbol {y})$. Integrating by parts twice leads to

(3.8)\begin{align} &\int_{\varOmega} \varphi \left[\overline{\boldsymbol{\nabla}_y \boldsymbol{\cdot} \left( \frac{1}{\rho_0}\,\boldsymbol{\nabla}_y H_{\varphi}\right) - \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{\rho_0^2 c_0^2}\,\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t} \right)}\right]\mathrm{d} \varOmega_y\nonumber\\ &\qquad+\int_{\partial\varOmega} \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{n}}{\rho_0 c_0^2} \left( \varphi\,\overline{\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t}} - \overline{H_{\varphi}}\,\frac{\mathrm{D}_0 \varphi}{\mathrm{D} t} \right) + \frac{1}{\rho_0} \left( \overline{H_{\varphi}}\,\frac{\partial \varphi}{\partial n} - \varphi\,\frac{\partial\overline{H_{\varphi}}}{\partial n} \right)\mathrm{d} \varGamma_y\nonumber\\ &\quad={-}\int_{\varOmega}\overline{H_{\varphi}}q\,\mathrm{d}{\varOmega_y} , \end{align}

where the integrals and the derivatives operate on the $\boldsymbol {y}$ coordinate.

Again, the objective is to tailor the Green's function $H_\varphi$ to obtain an explicit expression for $\varphi$ in terms of the sources $q$ and $v_s$. Requiring that the first integral in the above expression reduces to $\varphi (\boldsymbol {x})$ leads to the following adjoint equation for the Green's function:

(3.9)\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \frac{1}{\rho_0}\,\boldsymbol{\nabla} {H_{\varphi}}\right) - \rho_0\,{\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{\rho_0^2 c_0^2}\,\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t} \right)} = \delta(\boldsymbol{x}-\boldsymbol{y}) , \end{equation}

where the derivatives are taken with respect to $\boldsymbol {y}$. Comparing with the direct problem (3.1), it is clear that Pierce's equation is self-adjoint.

The boundary conditions for $H_\varphi$ are now selected to remove $\varphi$ and its derivatives from the boundary integral in (3.8). For the lined vibrating surface $\varGamma$, after using the fact that $\boldsymbol {u}_0\boldsymbol {\cdot }\boldsymbol {n}=0$ and introducing the Myers condition (3.3) for Pierce's equation, the contribution to the boundary integral becomes

(3.10)\begin{equation} I_\varGamma ={-}\int_\varGamma \overline{H_\varphi}\,\frac{\mathrm{D}_0\xi}{\mathrm{D} t} - \overline{H_\varphi}\xi\boldsymbol{n} \boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0] + \frac{\varphi}{\rho_0}\,\overline{\frac{\partial H_\varphi}{\partial n}}\,\mathrm{d}\varGamma . \end{equation}

Following the same reasoning as in § 2.2, this contribution can be written as

(3.11)\begin{equation} I_\varGamma = \int_{\varGamma}\frac{\rho_0 v_s}{\mathrm{i} \omega}\,\overline{\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t}} \,\mathrm{d} \varGamma , \end{equation}

provided that the tailored Green's function satisfies the following adjoint Myers condition for Pierce's equation:

(3.12)\begin{equation} \frac{\partial {H_{\varphi}}}{\partial n} = \rho_0^2\left\{\frac{\mathrm{D}_0\nu}{\mathrm{D} t} - \nu\boldsymbol{n}\boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0]\right\}, \quad\text{with}\ \nu =\frac{ \bar{{A}}}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0 }{\mathrm{D} t}\left(\frac{H_{\varphi}}{\rho_0} \right) . \end{equation}

On the far-field surface $\varGamma _{\infty }$, the solution $\varphi$ satisfies the radiation condition (2.16), and the adjoint radiation condition for $H_\varphi$ is readily identified as

(3.13)\begin{equation} \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{y}}{c_0^2 r}\,\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t} - \frac{\partial H_{\varphi}}{\partial r} ={-}\mathrm{i} k \bar{\alpha} H_{\varphi} , \end{equation}

which is the same as (2.30).

As a consequence, the potential $\varphi$ for an observer $\boldsymbol {x}$ can be written explicitly in terms of the sources $v_s$ and $q$ as

(3.14)\begin{equation} \varphi(\boldsymbol{x}) ={-} \int_{\varGamma} \frac{\rho_0}{\mathrm{i} \omega}\, \overline{\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t}}\,v_s \,\mathrm{d}{\varGamma_y} -\int_{\varOmega}\overline{H_{\varphi}} q \,\mathrm{d}{\varOmega_y} ={-}\left\langle\frac{\mathrm{D}_0 H_{\varphi}}{\mathrm{D} t},\frac{\rho_0 v_s}{\mathrm{i} \omega}\right\rangle_{\varGamma} - \left\langle H_{\varphi},q\right\rangle_{\varOmega} , \end{equation}

provided that the tailored Green's function $H_\varphi$ satisfies the adjoint problem defined by (3.9), (3.12) and (3.13).

All the comments in § 2.2.4 concerning the adjoint problem for Goldstein's equation also apply here. In particular, if one is interested in computing an acoustic quantity other than $\varphi$, then it is possible to modify the right-hand side of (3.9) accordingly. For instance, to obtain an explicit expression for the pressure field similar to (3.14), the required tailored adjoint Green's function $H_p$ satisfies

(3.15)\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \frac{1}{\rho_0}\,\boldsymbol{\nabla} {H_p}\right) - \rho_0\, {\frac{\mathrm{D}_0}{\mathrm{D} t} \left( \frac{1}{\rho_0^2 c_0^2}\,\frac{\mathrm{D}_0 H_p}{\mathrm{D} t} \right)} = \overline{\frac{\mathrm{D}_0}{\mathrm{D} t}}\,\delta(\boldsymbol{x}-\boldsymbol{y}) , \end{equation}

where the derivatives apply on the $\boldsymbol {y}$ coordinate. The tailored adjoint Green's function $H_p$ also satisfies the boundary conditions (3.12) and (3.13).

4.1. Flow reversal theorem

While the two propagation operators (2.1) and (3.1) considered in this paper are self-adjoint, the associated boundary conditions are not self-adjoint. As shown above, the Myers impedance condition (written in (2.5) and (3.3) for $\phi$ and $\varphi$, respectively) differs in the adjoint problems for the tailored Green's functions (see (2.27) and (3.12) for $\phi$ and $\varphi$, respectively). Similarly, the radiation conditions (2.16) and (2.15) differ between the direct and adjoint problems. This implies that non-reflecting conditions and buffer zones, such as perfectly matched layers, would have to be rewritten specifically for the adjoint problems where the radiation condition is anti-causal. Having to implement separate solvers for the direct and adjoint problems is a hindrance that can be partly avoided.

To that end, one has to reverse the mean flow direction (i.e. substitute $\boldsymbol {u}_0$ by $-\boldsymbol {u}_0$) and take the complex conjugate of the propagation equation and the boundary conditions. For Goldstein's equation, the governing equations (2.32), (2.27) and (2.16) for the tailored adjoint Green's function $G_\phi$ become

(4.1)\begin{equation} \left.\begin{aligned} & \rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t}\left(\frac{1}{ c_0^2 }\,\frac{\mathrm{D}_0\overline{G_{\phi}}}{\mathrm{D} t} \right) - \boldsymbol{\nabla} \boldsymbol{\cdot} \left(\rho_0\,\boldsymbol{\nabla} \overline{G_{\phi}} \right) = \delta(\boldsymbol x -\boldsymbol y) \quad \text{in }\varOmega,\\ & \frac{\partial {\overline{G_{\phi}}}}{\partial {n}} ={-}\frac{\mathrm{D}_0\bar{\eta} }{\mathrm{D} t} + \bar{\eta} \boldsymbol{n}\boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0]\quad\text{with } \bar{\eta} =\frac{\rho_0 {A}}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0\overline{G_{\phi}}}{\mathrm{D} t} \quad\text{on }\varGamma,\\ & \frac{\boldsymbol{u}_0\boldsymbol{\cdot}\boldsymbol{y}}{c_0^2 r}\,\frac{\mathrm{D}_0 \overline{G_{\phi}}}{\mathrm{D} t} - \frac{\partial \overline{G_{\phi}}}{\partial r} = \mathrm{i} k \alpha \overline{G_{\phi}}\quad\text{on } \varGamma_\infty. \end{aligned} \right\} \end{equation}

It is apparent that these equations for $\overline {G_\phi }$ are consistent with the direct problem for $\phi$ given in § 2.1. Hence existing solvers for the direct problem can be reused directly for the adjoint problem by solving for the complex conjugate of the adjoint Green's function and by reversing the mean flow direction. This so-called flow reversal theorem allows us to relate more easily the direct and adjoint problems; see Godin (Reference Godin1997), Möhring (Reference Möhring1978) and Eversman (Reference Eversman2001b).

For Pierce's equation, the governing equations (3.1), (3.6) and (2.16) for the tailored adjoint Green's function $H_\varphi$ become

(4.2)\begin{equation} \left. \begin{aligned} & \boldsymbol{\nabla} \boldsymbol{\cdot} \left(\frac{1}{\rho_0}\,\boldsymbol{\nabla} \overline{H_{\varphi}} \right) -\rho_0\,\frac{\mathrm{D}_0}{\mathrm{D} t}\left(\frac{1}{\rho_0^2 c_0^2}\,\frac{\mathrm{D}_0 \overline{H_{\varphi}}}{\mathrm{D} t} \right) = \delta(\boldsymbol{x}-\boldsymbol{y}) \quad\text{in }\varOmega, \\ & \frac{\partial {\overline{H_{\varphi}}}}{\partial {n}} ={-}\rho_0^2\left\{\frac{\mathrm{D}_0\bar{\nu} }{\mathrm{D} t} - \bar{\nu}\boldsymbol{n}\boldsymbol{\cdot}[(\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_0]\right\}\quad\text{with } \bar{\nu} =\frac{{A}}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0 }{\mathrm{D} t}\left(\frac{\overline{H_{\varphi}}}{\rho_0} \right) \quad\text{on }\varGamma,\\ & \frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{y}}{c_0^2 r}\,\frac{\mathrm{D}_0 \overline{H_{\varphi}}}{\mathrm{D} t} - \frac{\partial \overline{H_{\varphi}}}{\partial r} = \mathrm{i} k\alpha \overline{H_{\varphi}} \quad\text{on }\varGamma_\infty. \end{aligned} \right\} \end{equation}

It can be noted that the governing equation and the radiation condition are identical to the direct problem for $\varphi$; see (3.1) and (2.16). However, the impedance condition differs from that of the direct problem (3.6), except in the case of an incompressible base flow.

The variational form for the system (4.2) reads

(4.3)\begin{align} &\int_{\varOmega} \frac{1}{\rho_0}\left(\frac{1}{c_0^2}\,\frac{\mathrm{D}_0 \overline{H_{\varphi}}}{\mathrm{D} t}\,\overline{\frac{ \mathrm{D}_0 \psi}{\mathrm{D} t}}- \boldsymbol{\nabla} \overline{H_{\varphi}}\boldsymbol{\cdot} \boldsymbol{\nabla} \bar{\psi} \right) \mathrm{d} \varOmega\nonumber\\ &\quad -\int_{\partial\varOmega} \frac{\bar{\psi}}{\rho_0}\left(\frac{\boldsymbol{u}_0 \boldsymbol{\cdot} \boldsymbol{n}}{c_0^2}\,\frac{D_0 \overline{H_{\varphi}}}{\mathrm{D} t} - \frac{\partial \overline{H_{\varphi}}}{\partial n}\right) \mathrm{d} \varGamma = \int_{\varOmega}\bar{\psi} \delta(\boldsymbol{x}-\boldsymbol{y}) \,\mathrm{d}{\varOmega} . \end{align}

The contribution of the impedance condition to the boundary integral is

(4.4)\begin{equation} I_{\varGamma} = \int_{\varGamma}\frac{\rho_0{A}}{\mathrm{i} \omega }\, \frac{\mathrm{D}_0}{\mathrm{D} t}\left( \frac{\overline{H_{\varphi}}}{\rho_0}\right)\overline{\frac{\mathrm{D}_0 \psi}{\mathrm{D} t}}\,\mathrm{d} \varGamma . \end{equation}

4.2. High-order finite element model

We use a high-order continuous Galerkin finite element method (FEM) equipped with a basis of integrated Legendre polynomials to discretise the weak formulations of the direct and adjoint problems. In Bériot et al. (Reference Bériot, Prinn and Gabard2016), this approach was shown to provide substantial reductions in memory and CPU time when compared to conventional finite elements for acoustic applications. The benefits of a high-order FEM (e.g. low dispersion error, exponential convergence for smooth solutions) are also retained in the presence of background mean flows (Bériot, Gabard & Perrey-Debain Reference Bériot, Gabard and Perrey-Debain2013). In the numerical simulations, it is helpful to maintain an equivalent discretisation accuracy when varying the frequency or the Mach number. We resort to the a priori error indicator proposed in Bériot & Gabard (Reference Bériot and Gabard2019), which adjusts the order across the mesh so as to achieve a given, user-defined $L^2$-error target accuracy $E_T$. In practice, the edge orders are first determined based on a one-dimensional error indicator, which accounts for the local in-flow dispersion relation properties and possible edge curvature. In a second step, the element interior (directional) orders are assigned through a set of simple element-based dependent conformity rules. The orders are here defined to be in the range $p_{{FEM}} \in [1,10]$.

The governing equations of the direct and adjoint problems contain singular source terms, which require special attention. In practice, the point sources are enforced directly using a Dirac in the weak form. Error estimates have been derived in the literature for such elliptic problems with singular right-hand sides; see, for instance, Bertoluzza et al. (Reference Bertoluzza, Decoene, Lacouture and Martin2018). They indicate that the preponderant part of the error is located in the close vicinity around the source. A usual recommendation is therefore to use graded meshes, so as to confine the singularity errors in a more compact region. However, mesh refinements increase model complexity and may become unwieldy when several Dirac source terms are present. In Koppl & Wohlmuth (Reference Koppl and Wohlmuth2014), it is shown that the optimal convergence of high-order finite elements is recovered without mesh grading, if one excludes the one-ring neighbourhood elements to the point source from the error evaluation. In practice, one is usually not interested in getting a very accurate solution in the few elements directly surrounding the singularity. As a result, in this study, no mesh refinement is applied around the source when computing the tailored adjoint Green's functions.

4.3. Acoustic transfer vectors

Once discretised, the explicit expression (2.31) for the solution $\phi$ can be written as scalar products between complex-valued vectors:

(4.5)\begin{equation} \phi(\boldsymbol{x}) = \boldsymbol{a}_\phi^{\rm T}(\boldsymbol{x})\,\boldsymbol{v}_s + \boldsymbol{b}_\phi^{\rm T}(\boldsymbol{x})\,\boldsymbol{q} , \end{equation}

where $^{\rm T}$ denotes the Hermitian transpose. A similar expression can be written for $\varphi$ using (3.14). The vectors $\boldsymbol {v}_s$ and $\boldsymbol {q}$ contain the degrees of freedom representing the source terms $v_s$ and $q$ on the finite element mesh. The vectors $\boldsymbol {a}_\phi$ and $\boldsymbol {b}_\phi$ can be identified easily from (2.31) and computed using the solution for $G_\phi$. These vectors are sometimes referred to as acoustic transfer vectors as they act as transfer functions between the distributed sources $v_s$ and $q$ and the solution observed at a point $\boldsymbol {x}$ (Tournour et al. Reference Tournour, Cremers, Guisset, Augusztinovicz and MArki2000). These vectors are particularly useful when considering large numbers of different source distributions. For each individual source distribution, the new solution $\phi$ can be calculated very rapidly using (4.5) (a scalar product has linear complexity with respect to the number of degrees of freedom involved).

Additionally, when dealing with stochastic sources, acoustic transfer vectors can greatly simplify the computation of the acoustic field properties. For instance, for the volume source $q$, one can write

(4.6)\begin{equation} \mathbb{E}[\phi\bar{\phi}] = \boldsymbol{b}_\phi^{\rm T}\,\mathbb{E}[\boldsymbol{q}\boldsymbol{q}^{\rm T}]\,\boldsymbol{b}_\phi , \end{equation}

where $\mathbb {E}$ denotes the expected value. The correlation matrix $\mathbb {E}[\boldsymbol {q}\boldsymbol {q}^{\rm T}]$ can be either calculated or modelled, depending on the nature of the source mechanism. From this correlation matrix, the sound field can be computed directly using the acoustic transfer vector $\boldsymbol {b}_\phi$. This approach has been used extensively to predict aerodynamic noise generation by turbulent flows, starting from Tam & Auriault (Reference Tam and Auriault1999). See Spieser (Reference Spieser2020) for a comprehensive literature review on this approach.

It is important to note the presence of the gradient of $G_\phi$ in the definition of the vector $\boldsymbol {a}_\phi$. The gradient can be evaluated directly from the finite element approximation of $G_\phi$ obtained after solving (4.1). This can be an issue if a low-quality numerical solution is used for $G_\phi$ since the evaluation of $\boldsymbol {a}_\phi$ will be inaccurate.

An alternative approach is to integrate by parts the first integral in (2.31) to transfer the derivative onto $v_s$. This is similar to rewriting (2.7) into (2.10). The resulting expression is

(4.7)\begin{equation} \phi(\boldsymbol{x}) = \int_{\varGamma} \overline{G_{\phi}}\,\frac{\rho_0}{\mathrm{i} \omega}\,\frac{\mathrm{D}_0 v_s}{\mathrm{D} t}\,\mathrm{d}{\varGamma_y} -\int_{\varOmega}\overline{G_{\phi}}\,q \,\mathrm{d}{\varOmega_y} . \end{equation}

This avoids having to differentiate $G_\phi$, which is computed numerically. However, if the prescribed normal velocity $v_s$ is poorly represented, for instance with noisy data, then computing its derivative might introduce more numerical error. Depending on the application at hand, one has to choose which of these two approaches is best suited.

5.1. Cylinder in compressible flow

The proposed methodology is first applied to compute the sound radiated by a point source located in the vicinity of a cylinder immersed in a steady compressible flow. The computational domain consists of a square of side $8\ \mathrm {m}$, represented in figure 1(a). The cylinder is centred at $x=y=0$ and has radius $R=1\ \mathrm {m}$. Far from the cylinder, the mean flow is in the positive $x$ direction with Mach number $M_{\infty }=0.3$, density $\rho _{\infty }=1.2\ \mathrm {kg}\ \mathrm {m}^{-3}$, and speed of sound $c_{\infty }=340\ \mathrm {m} \mathrm {s}^{-1}$. The flow field around the cylinder is computed using a compressible homentropic potential flow solver (hence for an inviscid, non-heat-conducting gas). The flow is then interpolated onto the acoustic mesh. The local Mach number varies up to 0.662, as illustrated in figure 1(a). The cylinder is treated acoustically with uniform impedance $Z=\rho _{\infty }c_{\infty }(2-2\mathrm {i})$ and is vibrating with surface velocity $\boldsymbol {v}(x,y) = [\cos (x),0]^{\rm T}$. In this example, a point source $\delta (\boldsymbol {x}-\boldsymbol {x}_{s})$ is placed at $x_{s}=-1\ \mathrm {m}$ and $y_{s}=-2\ \mathrm {m}$, and the observer is located at $x=y=2.5\ \mathrm {m}$.

Figure 1. Computational domain: (a) mean Mach number, (b) mean density, and (c) mesh with FEM face order distribution for $\omega =10 c_0/R$.

The adjoint methods for Goldstein's and Pierce's equations are used to compute the acoustic field at the observer. Numerical simulations are performed over the normalised frequency range $1 \leqslant \omega R/c_{\infty } \leqslant 10$. The computational domain is discretised with 6-node triangular elements whose curved edges allow us to better represent the geometry of the cylinder. The size of the finite elements ranges from $h=0.05\ \mathrm {m}$ on the cylinder surface to $h=0.25\ \mathrm {m}$ at the outer boundary of the domain, corresponding respectively to 12.6 and 2.5 elements per acoustic wavelength at the highest frequency of interest. This leads to approximately 6700 elements in the domain. The FEM a priori error indicator automatically adjusts the order of the finite element basis to the user-defined accuracy, the frequency, the mean flow and the element size distribution, also accounting for the curvature of the mesh (Bériot & Gabard Reference Bériot and Gabard2019). The mesh is shown in figure 1(c) together with the distribution of polynomial orders obtained with target accuracy $E_T=1\,\%$. The face order $p_{{FEM}}$ varies from 2 to 5 at $\omega =10 c_{\infty }/R$. At the outer boundaries, perfectly matched layer (PML) regions made of ten layers are extruded automatically to absorb outgoing waves (Bériot & Modave Reference Bériot and Modave2021). They are stabilised using Lorentz transformation, following the work in Marchner et al. (Reference Marchner, Beriot, Antoine and Geuzaine2021).

Figure 2 presents the results obtained for the direct and adjoint problems at $\omega =10 c_{\infty }/R$ using Goldstein's equation. The real parts of the acoustic potential and pressure fields obtained for the direct problem are shown in figures 2(a) and 2(b), respectively. The sound field radiated from the source is clearly visible and is scattered by the cylinder, which also radiates sound due to its surface vibration. Also visible is the mean flow effect on the sound waves, with shorter waves propagating upstream, and longer waves propagating downstream. Figures 2(c) and 2(d) show the solutions $G_{\phi }$ and $G_{p}$ of the adjoint problem. In these solutions, waves are converging towards the observer location (these waves are anti-causal) and are scattered by the cylinder. Due to the flow reversal theorem, shorter wavelengths are observed on the right and longer wavelengths on the left, which is the opposite of the direct solutions.

Figure 2. Real parts of (a) the acoustic potential and (b) the pressure obtained with the direct problem for Goldstein's equation. Real parts of the numerical Green's functions (c) $G_{\phi }$ and (d) $G_{p}$ obtained from the corresponding adjoint problem, at $\omega =10 c_{\infty }/R$. Symbol $\times$ indicates source location, and $+$ indicates observer position.

The solutions of the direct problems based on Goldstein's and Pierce's equations are now compared to those obtained with the tailored adjoint Green's functions $G_\phi$, $G_p$ or $H_\varphi$, using (2.31), (2.33) or (3.14). Figure 3 shows the variables $\phi (\boldsymbol {x})$, $\varphi (\boldsymbol {x})$ and $p(\boldsymbol {x})$ computed at the observer location for a range of frequencies $1 \leqslant \omega R/c_{\infty } \leqslant 10$. For Goldstein's equation, the potential and pressure fields calculated with the adjoint approach and (2.31) and (2.33) match exactly the solutions of the direct problem. The same conclusion is obtained for Pierce's equation. These results demonstrate the applicability of the adjoint method to compressible non-parallel flows, in the presence of lined and vibrating surfaces.

Figure 3. Real parts of (a) the acoustic potential $\phi$, (b) the acoustic potential $\varphi$, and (c) the acoustic pressure computed at $x=y=2.5R$. Here, $\circ$ indicates Goldstein's direct problem, $\times$ indicates Pierce's direct problem, red solid line indicates Goldstein's adjoint problem, and dashed line indicates Pierce's adjoint problem.

For this test case, the assumptions of Goldstein's equation (homentropic and irrotational mean flow) are all valid, and this equation thus provides the reference solution. Note that the mean flow is not exactly isothermal since the mean temperature varies slightly around the cylinder (these variations do not exceed a few per cent). As a consequence, the mean density and sound speed also vary; see figure 1(b). In this case, Pierce's equation is approximate but still provides results in excellent agreement with Goldstein's equation (see figure 3c). It is only at low frequency ($\omega R/c_{\infty }\simeq 1$) that very small discrepancies can be observed between the two propagation equations. These discrepancies can be attributed to the small gradients of $\rho _0$ and $c_0$, and the high-frequency assumption used in Pierce's equation. In § 5.3, we will illustrate the differences between these wave equations when significant mean temperature gradients are present.

5.2. Sound source in a lined duct with flow

The second application involves a point source located in a two-dimensional straight duct with a uniform flow and a finite lined section modelled using the Myers impedance condition. The aim is to validate the formulation of the adjoint Myers condition for Goldstein's equation as well as the proposed formulations for the direct and adjoint Myers condition for Pierce's equation. These formulations, described in §§ 2 and 3, involve the neglect of contour terms around the lined surfaces, which could not be validated in the previous test case.

The straight duct, shown in figure 4(a), has height $h_{{d}}=0.2\ \mathrm {m}$ and length $r_{{d}}=0.75\ \mathrm {m}$. Inside the duct, the uniform flow is such that $\rho _0=1.225\ \mathrm {kg}\ \mathrm {m}^{-3}$, $c_0=340\ \mathrm {m}\ \mathrm {s}^{-1}$ and $M_0=u_0/c_0=0.5$. On the duct upper wall, uniform surface impedance $Z=\rho _0 c_0 (2-2\mathrm {i})$ is defined for $0.125< x<0.625$ m. The point source, placed at $x_s=0.25\ \mathrm {m}$ and $y_s=0.1\ \mathrm {m}$, is defined as $10^{-3} \delta (\boldsymbol x -\boldsymbol x_s )$. The acoustic potential and pressure fields are computed at the observer position $x=0.6\ \mathrm {m}$ and $y=0.1\ \mathrm {m}$.

Figure 4. (a) Computational domain with liner in grey. (b) Mesh with FEM face order distribution at $f=5$ kHz.

The finite element simulations are performed for frequency $f=5\ \mathrm {kHz}$ and target accuracy $E_T=1\,\%$. The mesh of triangular elements is shown in figure 4(b) together with the distribution of face orders. The size of the elements inside the duct varies from $h=0.0025\ \mathrm {m}$ at the wall to $h=0.025\ \mathrm {m}$ at the centre of the duct. The face order $p_{{FEM}}$ varies between 1 and 5. In order to absorb outgoing waves, PML regions made of 10 layers are introduced at the duct terminations.

Acoustic results obtained using Goldstein's equation are presented in figure 5. For the direct problem, the real part of the acoustic potential and pressure are given in figures 5(a) and 5(b), respectively. For the adjoint problem, the real part of the numerical Green's functions $G_{\phi }$ and $G_{p}$ are shown in figures 5(c) and 5(d). As for the previous test case, the Green's functions $G_{\phi }$ and $G_{p}$ exhibit waves converging towards the observer position, and the effect of the mean flow on the wavelength is reversed when compared to the direct problem. Using these tailored adjoint Green's functions together with (2.31) and (2.33) yields solutions for $\phi$ and $p$ at the observer location in figure 5 that are identical to the solutions of the direct problem within floating-point precision. The results obtained with Pierce's equation are not shown here for the sake of brevity since they are identical to those in figure 5. This is because, with a uniform mean flow, the wave operators in (2.1) and (3.1) are equivalent since dividing the former by $\rho _0^2$ yields the latter.

Figure 5. Real parts of (a) the acoustic potential and (b) the pressure field obtained for Goldstein's equation, and real parts of the corresponding numerical Green's functions (c) $G_{\phi }$ and (d) $G_{p}$ obtained from the adjoint problem at $f=5\ \mathrm {kHz}$. Symbol $\times$ indicates source location, and $+$ indicates observer location. The liner is indicated in grey.

To verify further the proposed models, they are compared to the semi-analytical mode-matching model from Gabard (Reference Gabard2010). The three solutions along the upper duct wall are compared in figure 6. Again, the solutions from Goldstein's and Pierce's equations are identical. They are also in close agreement with the mode-matching solution, even in the vicinity of the impedance discontinuities at $x=0.125$ and $0.625$, where the acoustic pressure varies rapidly. The velocity potentials $\phi$ and $\varphi$ vary continuously in these regions, but their slopes change across the impedance discontinuities. This behaviour is well resolved by the high-order finite element method. The pressure fields shown in figure 6 are obtained by computing $p=-\rho _0\,\mathrm {D}_0\phi /\mathrm {D} t$ or $p=\mathrm {D}_0\varphi /\mathrm {D} t$. It is the gradient $\boldsymbol {u}_0\boldsymbol {\cdot }\boldsymbol {\nabla }$ included in the material derivative that results in rapid variations of $p$ near impedance discontinuities. This behaviour is well-known and expected; see Gabard (Reference Gabard2010) and Luneville & Mercier (Reference Luneville and Mercier2014) for more details.

Figure 6. (a) Real part and (b) imaginary part of the pressure field at $f=5\ \mathrm {kHz}$ along the lined duct top wall ($y=h_{d}$). Here, blue dashed line indicates Pierce's equation, red solid line indicates Goldstein's equation, and $\circ$ indicates mode-matching results (reference).

These results confirm that the formulation of the adjoint Myers condition for Goldstein's equation proposed in § 2.2.2 is consistent with the standard formulation of Eversman (Reference Eversman2001a). This also confirms that the formulations of the direct and adjoint Myers impedance conditions for Pierce's equation proposed in § 3 are consistent with those for Goldstein's equation.

5.3. Propagation in sheared non-isothermal flows

In the results presented above, Pierce's and Goldstein's wave equations yield almost identical results. This is not always the case, in particular when the mean flow includes strong temperature gradients. To illustrate this, two additional test cases are considered, namely the sound radiated by a source in a hot jet, and the sound produced by a point source in a duct with a hot sheared flow. The solutions from the direct problems are also compared to reference solutions obtained with the linearised Euler equations solved using the same high-order finite element method (Hamiche et al. Reference Hamiche, Le Bras, Gabard and Beriot2019). The aim is to assess quantitatively the impact of the simplifications involved in Pierce's and Goldstein's operators when mean temperature gradients are present or when the mean flow is not homentropic.

Adjoint problems have also been computed and lead to the same results as the direct problems. For the sake of brevity, adjoint results have not been included in this subsection.

5.3.1. Free-field propagation

The case of a localised heat source radiating in a hot jet flow corresponds to the benchmark problem used in Spieser (Reference Spieser2020), except that two-dimensional calculations are performed instead of axisymmetric ones. The base flow is based on the fourth CAA workshop; see section ‘Radiation and refraction of sound waves through a two-dimensional shear layer’ of Dahl (Reference Dahl2004). The mean flow is parallel along the $x$ axis with a velocity profile defined as

(5.1)\begin{equation} u_{0}(y)= u_{j} \exp\left(-\frac{y^2}{2 \sigma^2}\right) , \end{equation}

where $\sigma\,{=}\,1.1018\ \mathrm {m}$, $u_{j}\,{=}\,M_{j} c_{j}$, $c_{j}\,{=}\,\sqrt {\gamma R T_{j}}$, $T_{j}\,{=}\,600\ \mathrm {K}$, $\gamma \,{=}\,1.4$ and $R\,{=}\,287\ \mathrm {m}^2\ \mathrm {s}^{-2}\ {\rm K}^{-1}$. The mean density is defined using the Crocco–Busemann law:

(5.2)\begin{equation} \frac{\rho_{j}}{\rho_0(y)} = \frac{T_{\infty}}{T_{j}} - \left(\frac{T_{\infty}}{T_{j}} -1\right)\frac{u_{0}(y)}{u_{j}} + \frac{\gamma-1}{2}\,M_{j}^2\, \frac{u_{0}(y)}{u_{j}} \left( 1-\frac{u_{0}(y)}{u_{j}}\right) , \end{equation}

where $T_{\infty }=300\ \mathrm {K}$ and $\rho _{j}=\sqrt {\gamma p_0/c_{j}^2}$, with $p_0=103\,330\ \mathrm {Pa}$. The jet flow Mach number is set to $M_{j} = 0.9$. Note that this mean flow is isentropic but not homentropic. Following Dahl (Reference Dahl2004), the sound source is defined as a source term in the energy equation of the LEE written for the variables $[\rho ',\rho _0 \boldsymbol u', p']^{\rm T}$. The source distribution is

(5.3)\begin{equation} q_p(x,y)= \sigma_s \sqrt{e} \exp\left( -\frac{x^2+y^2}{2 \sigma_s^2} \right) , \end{equation}

with $\sigma _s= 0.136\sigma$. Corresponding source terms are derived for Goldstein's and Pierce's equations.

The computational domain extends from $-5\sigma$ to $18\sigma$ in the $x$ direction, and from $-3\sigma$ to $9\sigma$ in the $y$ direction, with the jet centreline at $y=0$. At the external boundaries, a PML region of thickness $2\sigma$ and made of five layers is used to absorb outgoing waves. Simulations are performed for Strouhal numbers $St= 2\sigma f/u_{j}$ in the range $0.5 \leqslant St \leqslant 3$.

A triangular mesh is used for the finite element solution. The element size is $h=\sigma /30$ along the jet centreline to discretise the point source and the jet shear layers. It is increased progressively to reach $h=0.6\sigma$ at the external boundaries of the domain. For the FEM simulations, the target error for the a priori error estimator is set to $1\,\%$, resulting in FEM polynomial orders varying between 1 and 8.

The real part of the pressure field computed with Goldstein's, Pierce's and the LEE operators at $S_{t}=0.5$ is shown in figure 7. The real part of the pressure field obtained along the line $y=6\sigma$ is shown in figure 8(a). Some differences can be observed between the three models, especially downstream of the source. Pierce's equation follows more closely the LEE compared to Goldstein's equation. This is confirmed by computing a relative error with respect to the LEE solution as a function of frequency from the pressure field at $y=6\sigma$. Results are presented in figure 8(b). Results from Goldstein's equation are consistently less accurate than with Pierce's equation. Also significant in figure 8(b) is the decrease of the relative error as the frequency increases, which is consistent with the high-frequency analysis used by Pierce (Reference Pierce1990). This is explained by the fact that the assumptions needed for Goldstein's equation are not satisfied: the mean flow is neither homentropic nor irrotational. Pierce's equation does not need these assumptions, but relies instead on a high-frequency approximation.

Figure 7. Real part of the pressure field for the hot jet problem at $St=0.5$ obtained with (a) Goldstein's equation, (b) Pierce's equation, and (c) the LEE. Levels between $-10^{-5} \mathrm {Pa}$ and $10^{-5}\ \mathrm {Pa}$ from blue to yellow.

Figure 8. (a) Real part of the acoustic pressure at $y=6\sigma$ for $St=0.5$, and (b) corresponding relative error $\epsilon$ (with LEE as reference) for an increasing Strouhal number. Here, solid line indicates LEE, blue dotted line indicates Pierce's equation, and red dashed line indicates Goldstein's equation. Error levels for an isothermal jet with $T_0=600\ \mathrm {K}$ are also reported: $\times$ indicates Pierce's equation, and $\circ$ indicates Goldstein's equation.

Table 1 compares the computational costs of the different wave operators in terms of number of degrees of freedom (DoF), memory footprint and solving time. For a given mesh and accuracy, the use of the LEE leads to an increase by a factor 7 in the execution time and by a factor 15 in the memory requirements of the linear solver when compared to Pierce's or Goldstein's equations. This is consistent with the comparison for three-dimensional calculations reported by Hamiche et al. (Reference Hamiche, Le Bras, Gabard and Beriot2019).

Table 1. Comparison of the computational costs of the different wave operators as reported by the MUMPS linear solver (Amestoy et al. Reference Amestoy, Duff, L'Excellent and Koster2001), for the non-isothermal jet flow problem at $St=3$.

For completeness, figure 8(b) also includes results for the same case but with uniform temperature $T_0=600\ \mathrm {K}$. In this case, $\rho _0$ and $c_0$ are uniform, hence Goldstein's and Pierce's equations are equivalent, which is confirmed by figure 8(b). The accuracy of Goldstein's equation is improved significantly compared to the non-isothermal-jet case because the mean flow is now homentropic. The remaining error compared to the LEE is due to the jet flow, which contradicts the irrotational flow assumption used for Goldstein's equation.

5.3.2. In-duct propagation

In this test case, a point mass source is located in a duct containing a hot sheared flow. The definitions of the duct and the source are identical to those of the duct problem presented in § 5.2, except that no liner is considered. Inside the duct, a mean velocity profile is defined as

(5.4)\begin{equation} u_0(y) = u_{j} \left( 1-\left|1-\frac{2y}{h_{{d}}}\right|^n\right) , \end{equation}

with $n=3$, $T_{j}=600\ \mathrm {K}$, $\gamma =1.4$ and $R=287 \mathrm {m}^2\ \mathrm {s}^{-2}\ \mathrm {K}^{-1}$. The flow Mach number along the duct centreline is $M_{j} = 0.5$. The mean density is defined from the Crocco–Busemann law (5.2) with $T_{\infty }=300\ \mathrm {K}$ and $p_0=103\,330\ \mathrm {Pa}$.

The direct problems based on either Pierce's or Goldstein's operator are solved with the high-order FEM method, and the results are compared to the LEE solutions. All the numerical simulations are carried out with $f=5\ \mathrm {kHz}$ using a mesh similar to that used in § 5.2. The a priori accuracy target is $E_T=1\,\%$, resulting in polynomial orders $p_{{FEM}}$ varying between 1 and 4.

The acoustic pressure in the duct computed from the three linearised operators is shown in figure 9. While the overall pressure distributions are similar between the three solutions, the details of the solution obtained with Pierce's equation follow more closely the reference solution.

Figure 9. Real part of pressure in the duct with sheared flow at $f=5\ \mathrm {kHz}$, for (a) Goldstein's equation, (b) Pierce's equation, and (c) the LEE. Levels between $-7\ \mathrm {Pa}$ and $+7\ \mathrm {Pa}$ from blue to yellow.

The real and imaginary parts of the acoustic pressure along the duct top wall are shown in figure 10 for the three propagation operators. The best agreement with the LEE reference solution is clearly obtained with Pierce's equation. The relative error with respect to the LEE solution is $26.1\,\%$ using Pierce's equation and $52.3\,\%$ using Goldstein's equation. The variation of the relative error from the LEE solution is shown in figure 11 for a range of Helmholtz numbers $\omega h_{d}/c_{j}$ varying between 5 and 40. The same trend is observed where Pierce's equation is consistently more accurate than Goldstein's equation. These results indicate that Pierce's operator is more accurate when modelling acoustic propagation in regions with mean temperature gradients. Note that the peaks of error seen in figure 11 are due to the cut-off frequencies of the duct modes. Whenever the frequency of interest is close to the cut-off frequency of a duct mode, the models are very sensitive to any numerical error, and large discrepancies can be observed between different propagation models. This is discussed in more detail in § V.B. in Gabard (Reference Gabard2014).

Figure 10. (a) Real part and (b) imaginary part of the pressure field at $f=5\ \mathrm {kHz}$ along the duct top wall ($y=h_{d}$), Here, solid line indicates the LEE, blue dotted line indicates Pierce's equation, and red dashed line indicates Goldstein's equation.

Figure 11. Relative error with respect to the LEE solution for the acoustic pressure along the duct top wall ($y=h_{d}$). Here, blue dotted line indicates Pierce's equation, and red dashed line indicates Goldstein's equation.