Far-field approximations to the derivatives of Green’s function for the Ffowcs Williams and Hawkings equation

The surface correction to the quadrupole source term of the Ffowcs Williams and Hawkings integral in the frequency domain suffers from the computation of high-order derivatives of Green’s function. The far-field approximations to the derivatives of Green’s function have been used without derivation and verification in previous work. In this work, we provide the detailed derivations of the far-field approximations to the derivatives of Green’s function. The binomial expansions for the derivatives of Green’s function and the far-field condition are employed during the derivations to circumvent the difficulties in computing the high-order derivatives. The approximations to the derivatives of Green’s function are systemically verified by using the benchmarks two -dimensional convecting vortex and the co-rotating vortex pair. In addition, we provide the derivations of the approximations to the multiple integrals of Green’s function by using the far-field approximations to the derivatives.

where ρ = ρ − ρ o is the density fluctuation. ρ is the density, and ρ o is the ambient density at the far field. The far field pressure fluctuation p = p − p o is c 2 o ρ in the assumption that the acoustic wave is an isentropic process. p is the pressure, and p o is the ambient pressure at the far field. We use the superscript (·) to denote the fluctuation and subscript (·) o to denote the far-field variable. t is time.  (1) reduces to the classical form of FW-H equation [1,11] without boundary motions The solution to the FW-H equation can be expressed as surface integrals of the monopole source and dipole source and a volume integral of the quadrupole source. Hereinafter, we refer to these integrals as monopole source term, dipole source term, and quadrupole source term, respectively. The quadrupole source term is usually ignored under the assumption that the monopole source term and dipole source term dominate the far-field sound in low Mach number flows [12]. However, recent studies show that ignoring the quadrupole source term may result in significant spurious sound even at relatively low Mach numbers [13,14]. The spurious sound results from the vortical or entropy waves that pass across the integral surface [15]. Different methods have been proposed to eliminate or suppress the spurious sound [14,[16][17][18].
One of the most widely used methods for the elimination of spurious sound is the surface correction integral proposed by Wang et al. [19]. The idea of the surface correction integral is to account for the contribution of the quadrupole sources outside the integral domain because spurious sound is generated by the eddies crossing the boundary of the integral domain. With the frozen eddy assumption, the contribution of the quadrupole sources outside the integral domain can be modelled by using the quadrupole source flux on the boundary of the integral domain. Then, the surface correction integral can be constructed based on the quadrupole source flux accounting for the contributions of the quadrupole sources outside the integral domain. Spurious sound has been eliminated or dramatically reduced by the surface correction integral in the computation of the sound generated by flows over airfoils [19]. However, the surface correction integral requires that the eddies leave the boundary of the integration domain at a nearly constant speed. To improve the surface correction integral proposed by Wang et al. [19], different methods have been proposed based on accounting for the nonconstant convective velocity at the boundary of the integral domain [13,18]. Most of the surface correction integrals are developed within the time-domain framework.
Lockard and Casper [20] first proposed the surface correction integral for the FW-H equation within the frequency domain. The derivation of the surface correction integral for the frequency-domain method involves repeated integration by parts. Thus, the surface correction integral reported by Lockard and Casper [20] consists of a series of surface integrals as follows where f = 0 is the boundary of the volume integral of the quadrupole source in the FW-H equation (hereinafter referred to as the FW-H integral surface). i = √ −1 is the unit of the imaginary number. Without loss of generality, we set the convective velocity of the quadrupole sources along the o − x 1 direction and use U 1 to denote the convective velocity. G(x; y) is Green's function, where x and y are the observer and source positions, respectively. ω is the frequency. The surface correction integral for the frequency-domain method (Eq. (3)) works reasonably well in the test cases of a convecting vortex when the convective Mach number of the flow M a = U 1 /c o is lower than 0.5. The problems with the surface correction integral for the frequency-domain method (Eq. (3)) are that (1) the series of surface integrals on the right-hand side of Eq. (3) may diverge when the convective Mach number is higher than 0.5, and (2) the computation of the high-order derivatives of Green's function is quite complicated and nontrivial [20]. The divergence of Eq. (3) is caused by the convective velocity larger than the phase velocity of a perturbation propagating through convective flows. The far-field approximations to the derivatives of Green's function have also been utilized to fix the second problem of the surface correction integral. However, these far-field approximations to the derivatives of Green's function have been used without derivation and verification [21]. Detailed derivations and verification of the approximations to the derivatives of Green's function are given in this work. The derivations of the approximations are based on the binomial expansion of the derivatives of Green's function and the far-field condition. These approximations to the derivatives are verified in detail by computing the far-field acoustic pressure from the quadrupole source term generated by the two-dimensional convecting vortex and the co-rotating vortex pair. In addition, the approximations to the integrals of Green's function are derived by using the approximations to the derivatives of Green's function.
The remainder of the paper is organized as follows. We will give the derivations of the approximations to the derivatives of Green's function for 2D and 3D flows in Section 2.
The verifications of the approximations to the derivatives of Green's function are reported in Section 3. Finally, conclusions are drawn in Section 4.

Far-field approximations to the derivatives of Green's function
We investigate the FW-H equation in the form of a convective wave equation in the frequency domain in accordance with the work of Lockard et al. [10,20]. Green's function for the convective wave equation in the frequency domain is simply referred to as Green's function hereinafter when there is no ambiguity. We directly give the expressions of the far-field approximations to the derivatives of Green's function in subsection 2.1. Then, we provide the detailed derivations of the approximations in subsection 2.2.

Expressions of the far-field approximations to the derivatives of Green's function
The (l + 2)th-order derivative of Green's function with respect to the variables y q , y i and y j can be approximated at the far field as follows: where G is Green's function and ϕ is the phase function of Green's function.
For the two-dimensional (2D) flows, we have G ≈ G 2D and ϕ = ϕ 2D in the far field with where the observer and source locations are denoted by x and y, respectively. k = ω/c o is the wavenumber. β = √ 1 − M 2 is the Prantle-Glauert factor. M is the Mach number of the freestream flow. The distance R is computed by It is noted that the 2D Green's function G 2D given in Eq. (5) is the asymptotic Green's function for 2D flows. For the three-dimensional (3D) flows, we have G = G 3D and ϕ = ϕ 3D with where the distance R is computed by By using the approximations to the derivatives to Green's function (Eq. (4)), the surface correction integral (Eq. (3)) proposed by Lockard and Casper [20] can be approximated by The high-order derivatives of Green's function ∂ l−1 ∂y i ∂y j in Eq. (11). Therefore, the difficulties in computing the high-order derivatives of Green's function are circumvented.

Derivations of the approximations to the derivatives of Green's function at the far field
The approximations in Eq. (4) can be derived separately by setting q = 1, · · · , N D , where N D = 2 for two-dimensional flows and N D = 3 for three-dimensional flows. We report the detailed derivation of the approximations to the derivatives of Green's function when q = 1 in Eq. (4). The same derivations can be obtained following the same method for q = 2 and q = 3. The subscript q = 1 in Eq. (4) indicates that the lth-order derivative of ∂y i ∂y j is taken along the direction of the o − y 1 axis. When the subscript q = 1, Eq. (4) can be expressed as We re-write the left-hand-side of Eq. (12) as For the two-dimensional flows with G = G 2D , the lth-order derivative of Green's function with respect to the variable y 1 can be expressed with the binomial expansion as follows where C k 1 l is the binomial coefficient. We consider the far-field sound generated by the compact source flows. Without loss of generality, we set the origin of the reference frame near the compact source flows. For the far-field sound at the observation position x generated by the sources at the position y, we have |x| |y|. The k 1 th-order derivative to the distance R given by Eqs. (7) and (10) can be approximated by Similarly, we have The detailed derivations of Eqs. (15) and (16) are reported in Appendix B. By using Eq. (15), the k 1 th-order derivative of exp ϕ(x;y) can be approximated by The detailed derivation of Eq. (17) is reported in Appendix C. Thus, the binomial expansion of Eq. (5) can be expressed as Furthermore, the derivatives of R −1/2 terms on the right-hand side of Eq. (18) can be approximated according to Eq. (16). Thus, we obtain the far-field approximations to the (l + 2)th-order derivatives of Green's function for the two-dimensional flows with q = 1 as follows For three-dimensional flows with G = G 3D , the lth-order derivative of Green's function with respect to the variable y 1 can be expressed with the binomial expression as follows It is noted that Eq. (19) is not directly obtained from Eq. (18). A second-order derivative of Eq. (18) is first derived. Then, according to Eq. (16) and Appendix C, Eq. (19) can be obtained by ignoring the terms of high-order of R −1/2 and ϕ.
For the far-field sound at the observation position x generated by the sources at the position y, we have |x| |y|. The approximated k 1 th-order derivative to the distance R given by Eq. (15) is then employed to give the three-dimensional formulation of Eq. (17). Therefore, the binomial expansion as Eq. (20) can be expressed as Similar to the derivation of Eq. (15), the k 1 th-order derivative of R −1 can be approximated as follows Furthermore, the approximated derivatives of R −1 terms can be employed to simplify the right-hand side of Eq. (21). Thus, we obtain the far-field approximations to the (l + 2)th-order derivative of Green's function for the three-dimensional flows with q = 1 as follows The full derivations of the approximations given by Eq. (4) can be obtained by setting the subscript q = 2 for two-dimensional flows and q = 2, 3 for three-dimensional flows.

Results and discussion
We verify the proposed esimations by the benchmark cases of two-dimensional convecting vortex and co-rotating vortex pair. It is noted that both the benchmark cases reported in this section can be addressed by using the classical FW-H solver. Here, we have used different methods to solve the two benchmark cases to verify the formulations proposed in this manuscript.

Two-dimensional convecting vortex
We use the benchmark flow of a two-dimensional convecting vortex to verify the approximations to the derivatives of Green's function in this subsection. The time-varying velocity and pressure of the convecting vortex are given by: where p o is the pressure of the ambient medium at the far field. γ = 1.4 is the specific heat ratio of air. The quantities are normalized by using the ambient speed of sound c o , reference length L and density ρ of fluid. We assume that the vortex is convected along the o − y 1 axis with a velocity of U 1 . M a = U 1 /c o is the convecting Mach number, which is equal to the freestream Mach number. The parameters are given by a 0 = 1 and a 1 = 1/(2π), respectively. The parameters a 2 and r are computed by Figure 1 gives the schematic of the pressure field when the two-dimensional convecting vortex crosses the permeable FW-H integral surface. Here, we take the permeable FW-H integral surface as a square with an edge of 10 units in length, as shown by the dashed edges in Fig. 1. The contours in Fig. 1 show the distribution of the fluctuating pressure p − p o , which is normalized by p o .
The pressure fluctuation approaches zero monotonically away from the vortex centre, as shown in Fig. 2, where the convecting Mach number of the vortex is M a = 0.2. Thus, the sound pressure should approach zero at the far field [17]. However, the computation of the far-field sound pressure using the FW-H integral without the surface correction results in significant spurious sound when the vortex is moving across the permeable integral surface, as shown in Fig. 3, where the solid line shows the acoustic pressure computed by using the FW-H surface integral at the position x = (100, 0) that is completely spurious sound.
The surface correction integral given by Eq. (3) and its approximation at the far field given by Eq. (11) are shown by dashed lines and dash-dotted lines in Fig. 3, respectively. For convenient comparison, we plot the minus of the surface correction integral, i.e., −I Q , in Fig. 3. It is clear that the surface correction integrals computed with both the original derivatives of Green's function (Eq. (3)) and the approximated derivatives of Green's function (Eq. (11)) can accurately estimate the spurious sound generated by the vortex  3) and (11) in detail as follows.
We denote the first, second and third terms in Eq. (3) by and the first, second and third terms in Eq. (11) by where we use the subscript "appr" to indicate the approximated computation of the derivatives of Green's function in Eq. (11). The first terms in Eq. (3) and Eq. (11) are exactly the same since the approximation to the derivative of Green's function is not used. The differences between the second and third terms in Eq. (3) Fig. 4(a). The results indicate that the sound pressure p 2,appr approaches the sound pressure p 2 as the observer distance S increases. The difference between the sound pressure p 2,appr and p 2 is less than 5% of I Q (x) at the same position when S > 53, indicating that the right-hand side of Eq. (4) gives a reasonable approximation to its left-hand side when l = 1 when the observer distance S > 53. The variations in the sound pressure computed by using Eqs. (28) and (31) are shown in Fig. 4(b). Similar to that corresponding to the second term, the sound pressure p 3,appr approaches the sound pressure p 3 as the observer distance S increases. The difference between the sound pressure p 3,appr and p 3 is less than 5% of I Q (x) at the same position when S > 35, indicating that the right-hand side of Eq. (4) gives a reasonable approximation to its left-hand side when l = 2 when the observer distance S > 35. The results also indicate that the approximation on the left-hand side of Eq. (4) converges to its right-hand side faster with l = 2 than that with l = 1.  Fig. 5(a). Similar to Fig. 4(a), the sound pressure p 2,appr approaches the sound pressure p 2 as the observer distance S increases. The difference between the sound pressure p 2,appr and p 2 is less than 5% of I Q (x) at the same position when S > 71. The variations in the sound pressure computed by using Eqs. (28) and (31) are shown in Fig. 5(b). The sound pressure p 3,appr approaches the sound pressure p 3 as the observer distance S increases. The difference between the sound pressure p 3,appr and p 3 is less than 5% of I Q (x) at the same position when S > 56. The results also support the conclusion that the approximation on the left-hand side of Eq. (4) approaches its right-hand side more rapidly with l = 2 than with l = 1. Compared with those plotted in Fig. 4, the  Fig. 6(a). Similar to those plotted in Figs. 4(a) and 5(a), the sound pressure p 2,appr approaches the sound pressure p 2 as the observer distance S increases. The variations in the sound pressure computed by using Eqs. (28) and (31) are shown in Figure 6(b). The sound pressure p 3,appr approaches the sound pressure p 3 as the observer distance S increases. All of these results support the conclusion that the approximation on the left-hand side of Eq. (4) converges to its right-hand side faster with l = 2 than that with l = 1. Compared with the plots in Figs. 4 and 5, the results also indicate that the convergence of approximations on the lefthand side of Eq. (4) to its right-hand side becomes slow as the convecting Mach number increases.  Fig. 7(a). Similar to those plotted in Fig. 4(a), the sound pressure p 2,appr converges to the sound pressure p 2 as the observer distance S increases. The difference between the sound pressure p 2,appr and p 2 is less than 5% of I Q (x) at the same position when S > 42. The variations in the sound pressure computed by using Eqs. (28) and (31) are shown in Fig. 7(b). The sound pressure p 3,appr converges to the sound pressure p 3 as the observer distance S increases. The difference between the sound pressure p 3,appr and p 3 is less than 5% of I Q (x) at the same position when S > 39. For θ = 90 • , the variations in the sound pressure computed by  Fig. 8(a). The sound pressure p 2,appr converges to the sound pressure p 2 as the observer distance S increases. The difference between the sound pressure p 2,appr and p 2 is less than 5% of I Q (x) at the same position when S > 42. The variations in the sound pressure computed by using Eqs. (28) and (31) are shown in Fig. 8(b). The sound pressure p 3,appr converges to the sound pressure p 3 as the observer distance S increases. The difference between the sound pressure p 3,appr and p 3 is less than 5% of I Q (x) at the same position when S > 37. All of these results support the conclusion that the approximation on the left-hand side of Eq. (4) converges to its right-hand side faster with l = 2 than that with l = 1. Comparison of the results plotted in Figures 4, 7, and 8 indicates that the approximations on the left-hand side of Eq. (4) converge to its right-hand side and become faster as the directive angle θ increases.

Co-rotating vortex pair
In this subsection, we use the benchmark flow of a co-rotating vortex pair to verify the approximations to the derivatives of Green's function. The co-rotating vortex pair is a classical model for studying the jet engine noise generated from the interaction of coherent structures [22]. A co-rotating vortex pair can be modelled by two line vortices of equal strength separated by the distance 2d, as shown in Fig. 9 [23]. The sound generated by the co-rotating vortex pair can be computed by modifying Lighthill's acoustic analogy equation into the vortex sound equation as follows where the sound sources other than the divergence of the Lamb vector ω × v are neglected according to dimensional analysis [24]. The corresponding frequency-domain formulation can be expressed as Fig. 9 Schematics of a co-rotating vortex pair where the operator F defines the Fourier transform. Although we use Lighthill's equation here, the corresponding Green's function is consistent with Green's function of the FW-H equation without freestream flow. Using Gauss's theorem, the far-field acoustic pressure can be computed by For a co-rotating vortex pair, the Lamb vector ω × v can be expressed as where = 4πd 2 is the angular velocity of the vortex pair and ±y (ν) = ± (d cos( t), d sin( t)) represents the location of the vortex pair. By using the sifting property of the Dirac delta function, the frequency-domain far-field acoustic pressure generated by the vortex pair can be obtained by the acoustic pressure expressed as follows To verify the approximations to the derivatives of Green's function, we first conduct Taylor's expansion for Green's function, By using the approximations to the derivatives of the Green's function given by Eq. (4), we can approximate the expansion of Green's function as follows where G 0 and ∂ϕ 0 ∂y m are Green's function G and space derivatives of the phase function ∂ϕ ∂y m at the origin. By replacing the Green's function in Eq. (36) with the approximations in Eq.
(38), we can obtain the approximated acoustic pressure as follows It is noted that the first nonzero term in Eq. (39) is the quadrupole term, confirming that the quadrupole source term dominates the acoustic pressure generated by the co-rotating vortex pair [25]. Figure 10 compares the far-field acoustic pressure computed by using the approximated acoustic pressure formulation Eq. (39) with that computed by the acoustic pressure formulation Eq. (36). For the results reported in Fig. 10, the observer is located at ( 100 2π , 100 2π ), which is sufficiently far away from the centre of the corotating vortex pair. The fluid density is taken as 1kg/m 3 . The speed of sound is set as c = 100m/s, which is sufficiently large compared to the characteristic velocity of the vortex U = d. The distance between the two vortices is 2d = 0.2m, and the angular velocity is = 2π rad/s. The black solid line in Fig. 10 is the result computed by using the acoustic pressure expression of Eq. (36), while the red circle is computed by using the approximated acoustic pressure formula (Eq. (39) with n=2). It is observed that the approximated acoustic pressure formulation (Eq. (39)) gives a good approximation to the vortex acoustic pressure result (Eq. (36)), indicating that the approximations to the derivatives of Green's function given by Eq. (4) are correct.

Conclusion
Computation of the high-order derivatives of Green's function for the FW-H equation is required to eliminate the spurious sound associated with the quadrupole sources. The approximations to the derivatives of Green's function in the frequency domain have been used without derivation in previous work. This work provides the detailed derivations of the approximations to the derivatives of Green's function. The essential expressions and formulations associated with the binomial expression of the derivatives of Green's function and the far-field condition to obtain the approximations are provided in detail. The benchmark flows of the two-dimensional convecting vortex and co-rotating vortex pair are used to verify the approximations. The results show that the approximations can be accurate at different Mach numbers and observer directions as long as the distance is large enough. The derivations of the approximations to integrals of Green's function by using approximations to the derivatives are also reported in detail.

Appendix A: Approximations to the integrals of the Green's function
The approximations to the integrals of Green's function have been used by Zhou et al. [17] to fix the divergence problem of the frequency-domain surface correction integral [20]. However, the detailed derivations of the approximations to the integrals of Green's function have not been provided. In this Appendix, we show that the approximations to the integrals of Green's function can be derived based on the approximations to the derivatives of Green's function.
The approximations to the integrals of Green's function with respect to the variable y 1 at the far field used in the work of Zhou et al. [17] are given as follows where ∂ 2 G q,n (x; y) ∂y i ∂y j = I q,n ∂ 2 G(x; y) is the multiple integral of ∂ 2 G(x;y) ∂y i ∂y j with respect to the variable y q . For the ith integral (i < n) from the inner of Eq. (41), the independent variable is ξ i and the upper limit of the integral is ξ i+1 . When i = n, the independent variable is ξ n , and the upper limit of the integral is y q .
Equation (40) is obtained by proving the following equation We give the details of the derivation of Eq. (42) for the two-dimensional flows by using the mathematical induction method and the approximations of Eq. (4). The derivation of Eq. (42) for the three-dimensional flows can be obtained in a similar manner.
To use the mathematical induction method, we first prove that Eq. (42) is valid when l = 1 and q = 1. For l = 1 and q = 1, Eq. (41) reduces to After transforming the partial derivative with respect to y to x and using integration by parts, the right-hand side of Eq. (43) becomes Considering the limit of Green's function at infinity and transforming the partial derivative with respect to x back to y, Eq. (45) becomes proving that Eq. (42) is valid when l = 1 and q = 1.
According to the procedures of the mathematical induction method, we assume that Eq. (42) is valid when l = h and q = 1.
According to Eqs. (41) and (47), we have Using integration by parts, we reformulate Eq. (48) as follows  We use a homogeneous vortical flow with a uniform Lighthill stress tensor to verify the far-field approximations to the integrals of Green's function. We define the uniform Lighthill stress tensor as T ij = A cos(2πft), where A is 1/s 2 . The uniform vortical flow moves within a domain V of [ 0, 10] ×[ 0, 10]. The frequency f is taken as 1Hz. The speed of sound is taken as c 0 = 340m/s, and the density is taken as ρ = 1kg/m 3 . The farfield observer is located at (340000m, 340000m), which is sufficiently large compared to the wavelength of sound. The freestream Mach number is 0.3 along the o − y 1 axis. The contribution from the Lighthill stress tensor to the far-field sound can be computed by using the quadrupole source term in the FW-H integral as follows By using the far-field approximations to the integrals of Green's function (Eq. (40)), the far-field sound can be approximated by ∂y i ∂y j n q dy.
(54) Figure 11 shows the results calculated by using Eqs. (53) and (54), respectively. The results computed by using the surface integral (Eq. (54)) match well with the volume integral (Eq. (53)), indicating that the approximation to the integrals of Green's function is valid.
By defining λ 1 = g 1/2 1 r 1 , the 1st-to 3rd-order derivatives of the distance R can be reformulated as follows ∂R ∂y 1 = ∂R ∂λ 1 ∂λ 1 For the high-order derivatives of the distance R, we have We can estimate that g l ∼ O(1), since 0 ≤ β < 1, max(δ 1l ) = 1, and min(δ 1l ) = 0. For the computation of far-field sound, we approximately estimate the order of magnitude by using the assumption λ 1 ∼ O(R). We notice that this assumption is not true for the observer near the vertical direction. However, it gives a reasonable approximation to the order of magnitude for the observer at the most part of the region. The numerical result in Section 3 also shows that the above relations give good approximations to the far-field sound pressure. The detailed investigation to the effects of this assumption on the results is expected to be conducted in the future.
We denote R as DR m 1 λ 1 n 1 with D = 1, m 1 = 1 and n 1 = 0. Further, the derivative of is thus of the order O(R 1−k 1 ). Therefore, we can express the k 1 th-order derivative of R with respect to y 1 as It is noted that the first order derivative of R with respect to y 1 is ∂R in the most part of the region where λ 1 ∼ O(R). For Eq. (16), a similar formulation to the Eq. (63) can be obtained as follows The derivation of Eq. (61) can be obtained by a similar way to that of Eq. (59) with the use of chain rule. It is noted that R −1/2 can be written as DR m 1 λ 1 n 1 with D = 1, The detailed derivation of Eq. (17) is reported in this Appendix. We expand the k 1 thorder derivative of exp ϕ as follows where E 1 , E 2 , E 3 , E 4 ... are the coefficients corresponding to the number of derivative's order k 1 . We start from the derivation from the Eq. (17) in the two-dimensional space. According to Eq. (6), the exponent ϕ of the 2D Green's function for the FW-H equation is as follows where i = √ −1 is the unit of the imaginary number. k = ω/c o is the acoustic wavenumber where ω and c o are the angular frequency and the speed of sound, respectively. x and y are the coordinate of observer and source locations, respectively. β = √ 1 − M 2 is the Prantle-Glauert factor. M is the Mach number of the freestream flow.
The derivative of the exponent ϕ with respect to y 1 is The k 1 th-order (k 1 ≥ 2) derivative of ϕ with respect to y 1 is Equation (17) in Section 2 can be obtained by combining Eqs. (67) and (69). The derivation of Eq. (17) for in three-dimensional space can be obtained in a similar way.