Phase-preserving narrow- and wide-angle parabolic equations for sound propagation in moving media

Ostashev, Vladimir E.; Colas, Jules; Dragna, Didier; Wilson, D. Keith

doi:10.1121/10.0024460

Parabolic equations are among the most popular numerical techniques in many fields of physics. This article considers extra-wide-angle parabolic equations, wide-angle parabolic equations, and narrow-angle parabolic equations (EWAPEs, WAPEs, and NAPEs, respectively) for sound propagation in moving inhomogeneous media with arbitrarily large variations in the sound speed and Mach number of the (subsonic) wind speed. Within their ranges of applicability, these parabolic equations exactly describe the phase of the sound waves and are, thus, termed the phase-preserving EWAPE, WAPE, and NAPE. Although variations in the sound speed and Mach number are often relatively small, omitting the second-order terms pertinent to these quantities can result in large cumulative phase errors for long propagation ranges. Therefore, the phase-preserving EWAPE, WAPE, and NAPE can be preferable in applications. Numerical implementation of the latter two equations can be performed with minimal modifications to existing codes and is computationally efficient. Numerical results demonstrate that the phase-preserving WAPE and NAPE provide more accurate results than the WAPE and NAPE based on the effective sound speed approximation.

I. INTRODUCTION

The narrow-angle parabolic equation (NAPE) and wide-angle parabolic equation (WAPE) are among the most popular computational techniques in atmospheric acoustics,^1–5 ocean acoustics,^6–8 nonlinear acoustics,^9,10 biomedical acoustics,¹¹ and electromagnetic propagation.¹² The NAPE and WAPE are well suited to relatively small computers, large domains, and high frequencies and can handle many phenomena such as atmospheric and ocean stratification and refraction; scattering by turbulence, internal waves, and other inhomogeneities; ground impedance and ocean bottom interactions; and propagation over slowly varying terrain and ocean bathymetry.

In atmospheric acoustics, the NAPE and WAPE should account for the effect of the wind velocity on the sound field. The corresponding equations have been formulated in a number of references, e.g., Ref. 2 and references therein. The WAPEs used in atmospheric acoustics are rather involved because they include the derivatives of the sound, density, and wind velocity, e.g., Refs. 13 and 14. These NAPEs and WAPEs are usually formulated assuming that the variations in the sound speed are small and the wind-speed Mach numbers are low.

A recent article¹⁵ suggests considering EWAPEs (extra-wide-angle parabolic equations or one-way equations), NAPEs, and WAPEs in the high-frequency approximation, in which the derivatives of the ambient quantities can be omitted. (Also, see Sec. 11.2 in Ref. 2.) The resulting PEs become much simpler than the previous formulations and still can be used in many applications. Specifically, Ref. 15 considers, in detail, the EWAPE given by Eq. (39) and the NAPE and WAPE derived from that equation. The WAPE from Ref. 15 has already been used in the literature, e.g., Refs. 16–19.

Reference 15 also suggests the other EWAPE, which is given by Eq. (B1), but considers the latter equation very briefly. One of the goals of the present article is to study this EWAPE in detail and determine when it is preferable compared to the EWAPE given by Eq. (39). To this end, starting with the EWAPE given by Eq. (B1) in Ref. 15, the NAPE and WAPE are derived. The newly formulated NAPE is remarkably simple and valid for arbitrary variations in the sound speed and Mach numbers. The new NAPE and the WAPE in the Padé (1,1) approximation can be implemented numerically with minimal modifications to existing codes and are computationally efficient.

In the literature, EWAPEs, NAPEs, and WAPEs are derived by omitting some terms in the Helmholtz-type equations. The ranges of applicability of the resulting equations are often difficult to assess. Hence, another goal of the present article is to study the ranges of applicability of these equations, which is examined by two approaches. In the first approach, starting with these equations, the dispersion relations (which are equivalent to the eikonal equations) are derived and compared to geometrical acoustics. When these dispersion relations agree with geometrical acoustics, the corresponding parabolic equations exactly describe the phase of a sound wave. The present article formulates the EWAPE, NAPE, and WAPE exactly describing the phase of a sound wave within the ranges of their applicability for arbitrary variations in the sound speed and Mach numbers. Therefore, these equations are termed the phase-preserving parabolic equations. On the other hand, many NAPEs and WAPEs from the literature preserve the phase of a sound wave only for low Mach numbers and/or small variations in the sound speed when the second-order terms can be omitted. Although the variations in the sound speed and Mach numbers are often relatively small, omitting the second-order terms might result in large cumulative phase errors for long ranges. For example, Ref. 15 argues that the second-order small terms cannot be omitted for the wind speed greater than or about 15 m/s; in this case, equations for arbitrary Mach numbers should be used. Thus, in many applications (such as sound propagation in the near-ground atmosphere and infrasound propagation in the upper atmosphere), it is desirable to use the phase-preserving EWAPE, NAPE, and WAPE.

The second approach for studying the ranges of applicability of parabolic equations is to compare them with the exact equation for sound propagation in stratified moving media. This article shows that the phase-preserving EWAPE, NAPE, and WAPE correctly describe sound propagation in stratified moving media within the ranges of their applicability.

This article is organized as follows. Section II briefly outlines results from Ref. 15 and presents new results such as the phase-preserving NAPE. In Sec. III, the ranges of applicability of EWAPEs, NAPEs, and WAPEs are studied by comparing them to the geometrical acoustics equations. In Sec. IV, these parabolic equations are compared to the exact equation in a stratified moving medium. Section V provides numerical algorithms for the phase-preserving NAPE and WAPE. Numerical results are presented in Sec. VI. The results are summarized in Sec. VII.

II. PARABOLIC EQUATIONS

This section considers parabolic equations for sound propagation in moving media.

A. EWAPE

First, we briefly overview some results from Ref. 15 and provide additional analysis of these results. Reference 15 expresses the sound pressure p of a monochromatic sound wave in a moving medium in terms of the auxiliary function

ϕ

such that

p (R) = \sqrt{\frac{ϱ}{ϱ_{0}}} (1 + \frac{i}{ω} v \cdot \nabla) ϕ (R) .

(1)

Here,

R = (x, y, z)

is the vector of Cartesian coordinates,

v (R)

is the medium velocity,

ϱ (R)

and

ϱ_{0}

are the density and its reference value, respectively, and ω is the angular frequency of the sound wave. The auxiliary function satisfies the convective Helmholtz equation,

[Δ - \frac{{(i ω - v \cdot \nabla)}^{2}}{c^{2}}] ϕ (R) = 0,

(2)

where

c (R)

is the sound speed. Equations and are derived by applying the high-frequency approximation (for which the sound wavelength is smaller than the characteristic scale of medium inhomogeneities) to the linearized equations of fluid dynamics. This approach was originally suggested in Ref. 20. Additionally, Eq. is derived assuming that

\nabla \cdot v = 0

⁠. This equality is valid in a stratified moving medium. It is also approximately valid if the deviation from a stratified medium is relatively small or the derivatives of the medium velocity can be omitted.

Starting with Eq. , the following one-way equation is obtained [Eq. (38) in Ref. 15]:

(\frac{\partial}{\partial x} + i k_{c} \hat{τ} - i σ k_{0} γ_{x} \sqrt{1 + {\hat{μ}}_{σ} + {\hat{χ}}_{σ}}) ϕ = 0.

(3)

Equation describes sound propagation in the positive direction of the x axis and, following Refs. 15 and 21, is termed the EWAPE (or the parabolic equation for extra-wide angles). In Eq. , σ is a factor specified below and

k_{0} = \frac{ω}{c_{0}}, k_{c} = \frac{ω}{c}, M_{x} = \frac{v_{x}}{c}, γ_{x} = \frac{1}{\sqrt{1 - M_{x}^{2}}}, {\hat{m}}_{⊥} = \frac{i}{ω} v_{⊥} \cdot \nabla_{⊥}, \hat{τ} = M_{x} γ_{x}^{2} (1 + {\hat{m}}_{⊥}) .

(4)

Here, k₀ and c₀ are the reference wavenumber and sound speed, respectively; k_c is the wavenumber in a motionless medium; v_x and

v_{⊥} = (v_{y}, v_{z})

are the velocity components in the direction of the x axis and the transverse plane, respectively; M_x and γ_x are the Mach number and Lorentz factor, respectively, pertinent to v_x;

\nabla_{⊥} = (\partial / \partial y, \partial / \partial z)

⁠; and the operators

{\hat{m}}_{⊥}

and

\hat{τ}

act on the transverse coordinates

r = (y, z)

⁠. Furthermore, in Eq. ,

{\hat{μ}}_{σ} = \frac{1}{σ^{2} k_{0}^{2}} Δ_{⊥}, {\hat{χ}}_{σ} = \frac{k_{c}^{2} γ_{x}^{2}}{σ^{2} k_{0}^{2}} {(1 + {\hat{m}}_{⊥})}^{2} - 1,

(5)

where

Δ_{⊥} = (\partial^{2} / \partial y^{2}, \partial^{2} / \partial z^{2})

⁠. Note that the factor σ appears in the notations given by Eq. but not in Eq. .

Equation is valid if the angle θ between the direction of sound propagation and the x axis is less than 90°. In deriving this equation, the derivatives of c and v are omitted; these derivatives are beyond the measurement accuracy in many applications. It is important to emphasize that Eq. is valid for arbitrary subsonic Mach numbers,

M = v / c

⁠, and arbitrary variations in the sound speed (provided that the derivatives of c and v can be omitted), which are characterized by the function ε:

ε = \frac{c_{0}^{2}}{c^{2}} - 1 = \frac{k_{c}^{2}}{k_{0}^{2}} - 1.

(6)

The function ε is the deviation of the refractive index squared in a motionless medium from one. In contrast, many parabolic equations from the literature are obtained assuming that ε and/or M are small parameters, keeping terms of order ε and M, and omitting higher-order terms such as

ε^{2}

and M². However, omitting the second- and higher-order terms may result in significant phase errors for relatively long propagation ranges. For example, Ref. 15 argues that for the wind speed greater than 15 m/s (⁠

M \sim 0.05

⁠), the second- and higher-order terms in M should not be neglected such that parabolic equations valid for arbitrary Mach numbers should be used.

The EWAPE given by Eq. does not actually depend on the factor σ. (This can be readily shown by moving σ premultiplying the square root into this square root.) However, when a NAPE or WAPE is obtained from Eq. by approximating the square-root operator in this equation, it does depend on the choice of σ. We desire that the resulting NAPE and WAPE be valid for arbitrary variations in the sound speed and Mach numbers, respectively, i.e., arbitrary ε and M. First, consider the NAPE, which can be obtained from Eq. by approximating the square-root operator with the first two terms of the Taylor series. This approximation can be used only if the terms containing the operators

{\hat{μ}}_{σ}

and

{\hat{χ}}_{σ}

in Eq. are small compared to one. To estimate these operators, we use a zeroth-order approximation for the auxiliary function:

ϕ \sim exp [i k_{0} (x cos θ + e \cdot r sin θ)]

⁠, where e is a unit vector characterizing the propagation direction of the sound wave in the (y,z)-plane. For relatively small angles θ, for which the NAPE is applicable, we have

| \nabla_{⊥} ϕ | \sim k_{0} θ ϕ

and

Δ_{⊥} ϕ \sim k_{0}^{2} θ^{2} ϕ

⁠. Then, in Eq. , the operators

{\hat{m}}_{⊥}, {\hat{μ}}_{σ}

⁠, and

{\hat{χ}}_{σ}

can be estimated as

{\hat{m}}_{⊥} \sim M θ, {\hat{μ}}_{σ} \sim \frac{θ^{2}}{σ^{2}}, {\hat{χ}}_{σ} \sim \frac{k_{c}^{2} γ_{x}^{2}}{σ^{2} k_{0}^{2}} {(1 + M θ)}^{2} - 1.

(7)

It follows from these relationships that for small θ, the terms containing the operator

{\hat{μ}}_{σ}

are small compared to one. However, the terms containing the operator

{\hat{χ}}_{σ}

are small only if

σ^{2} k_{0}^{2} = k_{c}^{2} γ_{x}^{2}

or

σ = \frac{k_{c} γ_{x}}{k_{0}} .

(8)

This is the desired value of σ for which the square-root operator in Eq. can be accurately approximated with the first two terms of the Taylor series for arbitrary ε and M.

Next, consider the WAPE obtainable from Eq. (3) by approximating the square-root operator with the Padé series expansion. The Padé series converges for all values of the operators ${\hat{μ}}_{σ}$ and ${\hat{χ}}_{σ}$ and, hence, for all values of σ.²² However, if ${\hat{μ}}_{σ}$ and ${\hat{χ}}_{σ}$ are small compared to one, fewer terms in this series are needed for a good approximation of the square root than in the case of large ${\hat{μ}}_{σ}$ and ${\hat{χ}}_{σ}$ ⁠. (This result is illustrated in Table I of Ref. 23.) This, again, dictates the choice of σ given by Eq. (8), at least for small θ.

TABLE I.

Comparison between the NAPEs and WAPEs. The first two columns specify the parabolic equation and the function which satisfies this equation. The coefficients $h_{n, m}$ pertinent to these equations are in the last four columns.

Equation	Function	$h_{1, 0}$	$h_{1, 2}$	$h_{2, 0}$	$h_{2, 2}$
NAPE	$\hat{ϕ}$	1	0	$n_{eff} - 1$	$\frac{1}{2 n}$
WAPE	$\hat{ϕ}$	1	$\frac{1}{4 n^{2} γ_{x}^{2}}$	$n_{eff} - 1$	$\frac{1}{2 n} + \frac{n_{eff} - 1}{4 n^{2} γ_{x}^{2}}$
WAPE_M	$\hat{ϕ}$	$1 + \frac{n^{2} - 1}{4}$	$\frac{1}{4 γ_{x}^{2}}$	$\frac{(n^{2} - 1) γ_{x}^{2}}{2} - (1 + \frac{n^{2} - 1}{4}) τ$	$\frac{1}{2} - \frac{τ}{4 γ_{x}^{2}}$
NAPE_eff	$\hat{p}$	1	0	$\frac{n_{eff}^{2} - 1}{2}$	$\frac{1}{2}$
WAPE_eff	$\hat{p}$	$1 + \frac{n_{eff}^{2} - 1}{4}$	$\frac{1}{4}$	$\frac{n_{eff}^{2} - 1}{2}$	$\frac{1}{2}$

Thus, the factor σ given by Eq. is preferable for approximating the square-root operator in the EWAPE [Eq. ]. Substituting Eq. into Eqs. and , we obtain

(\frac{\partial}{\partial x} + i k_{c} \hat{τ} - i k_{c} γ_{x}^{2} \sqrt{1 + \hat{μ} + \hat{χ}}) ϕ = 0,

(9)

\hat{μ} = \frac{1}{k_{c}^{2} γ_{x}^{2}} Δ_{⊥}, \hat{χ} = 2 {\hat{m}}_{⊥} + {\hat{m}}_{⊥}^{2} .

(10)

Equation is the desired EWAPE with the operators

\hat{μ}

and

\hat{χ}

given by Eq. . In these operators, the subscript σ is omitted because the factor σ is specified with Eq. . With σ given by Eq. and for small θ, the operators

\hat{μ}

and

\hat{χ}

are estimated as

\hat{μ} \sim \frac{k_{0}^{2} θ^{2}}{k_{c}^{2} γ_{x}^{2}}, \hat{χ} \sim M θ + {(M θ)}^{2} .

(11)

The terms containing these operators are small compared to one for small θ. This leads to the NAPE that is valid for arbitrary ε and M and requires fewer terms in the Padé series expansion for the WAPE.

Equations (9) and (10) coincide with Eqs. (B1) and (B2) in Ref. 15. However, Ref. 15 considers the EWAPE given by Eq. (9) only briefly in Appendix B. One of the main goals of this article is to study in detail this EWAPE and the NAPE and WAPE derived from this equation.

In the literature, the factor σ differs from that given by Eq. (8). In a motionless medium, σ is often set to one, e.g., Ref. 7. In this case, ${\hat{χ}}_{σ} = ε$ and the Taylor series expansion of the square-root operator in Eq. (3) is valid only for small ε. In a moving medium, the bulk of Ref. 15 considers the EWAPE, NAPE, and WAPE pertinent to $σ = γ_{x}$ ⁠. With this choice of σ, the NAPE is valid only for small ε, whereas the WAPE requires more terms in the Padé series.

Thus, the EWAPE given by Eq. (9) might be preferable compared to other equations used in the literature because the terms containing the operators $\hat{μ}$ and $\hat{χ}$ are small compared to one for arbitrary ε and M (for small θ). This feature of Eq. (9) is achieved by the factor $k_{c} γ_{x}^{2}$ premultiplying the square-root operator, whereas other EWAPEs do not have this factor. The EWAPE given by Eq. (9) does not contain the reference wavenumber k₀. It is worthwhile to point out that Eqs. (9) and (10) are simpler (i.e., have fewer terms) than the corresponding Eqs. (39) and (40) in Ref. 15, which pertain to $σ = γ_{x}$ ⁠. As expected, these equations coincide if c is constant.

In the literature, EWAPEs, NAPEs, and WAPEs are derived by omitting some terms in starting equations such as Eq. (2). The ranges of applicability of the resulting equations are often difficult to assess. Another goal of this article is to study the ranges of applicability of these equations in detail. To this end, two approaches are used. First, starting with the EWAPE, NAPE, and WAPE, the dispersion relations (which are equivalent to the eikonal equations) are derived and compared to geometrical acoustics. Section III shows that within the ranges of applicability of the EWAPEs given by Eqs. (3) and (9) (⁠ $θ < 90 °$ ⁠), the dispersion relations following from these equations coincide with geometrical acoustics. Therefore, Eqs. (3) and (9) are termed the phase-preserving EWAPEs.

In the second approach, EWAPEs, NAPEs, and WAPEs are compared with the exact equation for sound propagation in stratified moving media. Section IV shows that for $θ < 90 °$ ⁠, the EWAPEs given by Eqs. (3) and (9) coincide with this exact equation.

Equations presented in this section are the starting point of the subsequent analysis.

B. NAPE

To derive the NAPE from Eq. , the operator

\hat{χ} = 2 {\hat{m}}_{⊥} + {\hat{m}}_{⊥}^{2}

is substituted into this equation. Then, we take into account that in the square-root operator in the resulting equation,

{\hat{m}}_{⊥} \sim M θ

while

\hat{μ} \sim k_{0}^{2} θ^{2} / (k_{c}^{2} γ_{x}^{2})

⁠; see Eqs. and . Approximating the square-root operator with the Taylor series and keeping terms of order θ and

θ^{2}

yields

\sqrt{1 + \hat{μ} + \hat{χ}} = \sqrt{1 + \hat{μ} + 2 {\hat{m}}_{⊥} + {\hat{m}}_{⊥}^{2}} \approx 1 + \hat{μ} / 2 + {\hat{m}}_{⊥} .

(12)

Note that in this formula, terms of order

{\hat{m}}_{⊥}^{2}

cancel out. Substituting Eq. into Eq. , we obtain

(\frac{\partial}{\partial x} - \frac{i}{2 k_{c}} Δ_{⊥} - \frac{i ω}{c_{eff}} + \frac{1}{c_{eff}} v_{⊥} \cdot \nabla_{⊥}) ϕ = 0.

(13)

Here,

c_{eff} = c + v_{x}

is the effective sound speed, which appears in this equation as a convenient notation rather than an approximation. Equation , which is a new result, is the desired NAPE that is valid for the arbitrary variations in the sound speed and arbitrary Mach numbers, i.e., arbitrary ε and M. Equation does not contain the reference wavenumber k₀. After the auxiliary function

ϕ

is calculated with Eq. , the sound pressure is determined with Eq. . In a motionless medium, Eq. coincides with Eq. (3.38) from Ref. 6 and Eq. (4.17) from Ref. 24 if in these equations, the derivatives of the sound speed are set to zero.

The NAPE [Eq. (13)] is valid for relatively small propagation angles θ. By studying the phase errors of a sound wave in motionless and moving media, Refs. 7 and 15 conclude that NAPEs are valid for the propagation angles $θ ≲ 20 °$ ⁠. However, this conclusion depends on an acceptable phase error and the problem considered. For near-ground sound propagation, there is generally a consensus^2,13 that a NAPE [including Eq. (13)] can be used if $θ ≲ 20 °$ ⁠.

Section III shows that the dispersion relation derived from Eq. (13) coincides with geometrical acoustics within the ranges of applicability of Eq. (13); therefore, this equation is termed the phase-preserving NAPE. Section IV shows that within the ranges of its applicability, Eq. (13) coincides with the exact equation in a stratified moving medium. Given the wide range of applicability of the phase-preserving NAPE [Eq. (13)], it is remarkably simple. Numerical implementation of Eq. (13) is provided in Sec. V A.

Other NAPEs appearing in the literature are valid either for low Mach numbers and/or small variations in the sound speed. They preserve the phase of a sound wave only approximately and do not exactly coincide with the exact equation in a stratified moving medium. For example, the NAPE given by Eq. (46) in Ref. 15 is valid only for $| ε | ≪ 1$ ⁠; that equation is also much more involved than the phase-preserving NAPE [Eq. (13)].

The NAPE in a moving medium is usually formulated for the sound pressure p rather than for the auxiliary function

ϕ

⁠. It is often written in the form

[\frac{\partial}{\partial x} - \frac{i}{2 k_{0}} Δ_{⊥} - \frac{i k_{0}}{2} (1 + \frac{c_{0}^{2}}{c_{eff}^{2}}) + \frac{1}{c_{0}} v_{⊥} \cdot \nabla_{⊥}] p = 0.

(14)

[Equation coincides with Eq. (2.111) in Ref. 2 if p is replaced with

\hat{p} exp (i k_{0} x)

⁠.] Equation describes sound propagation only for small variations in the sound speed ε and low Mach numbers M. With the last term omitted, Eq. can be obtained from the NAPE in a motionless medium using the effective sound speed approximation; that is, by replacing c with

c_{eff}

⁠.¹ Because Eq. is valid only for small ε and M the deviation of

c_{eff}

from c₀ is also relatively small, i.e.,

| c_{eff} - c_{0} | / c_{0} ≪ 1

⁠. Therefore, keeping

c_{eff}^{2}

in Eq. is inconsistent with its ranges of applicability. It can be shown that within these ranges,

1 + c_{0}^{2} / c_{eff}^{2}

in Eq. can be replaced with

2 c_{0} / c_{eff}

⁠. In this case, the third terms in Eqs. and coincide.

C. WAPE

WAPEs are usually derived by approximating the square-root operator in Eq. with the Padé (m,m) series expansion such that

\sqrt{1 + \hat{μ} + \hat{χ}} = 1 + \sum_{j = 1}^{m} \frac{a_{j, m} (\hat{μ} + \hat{χ})}{1 + b_{j, m} (\hat{μ} + \hat{χ})} .

(15)

In this series, the coefficients

a_{j, m}

and

b_{j, m}

are given by²⁵

a_{j, m} = \frac{2}{2 m + 1} {sin}^{2} \frac{j π}{2 m + 1}, b_{j, m} = {cos}^{2} \frac{j π}{2 m + 1} .

(16)

Substituting Eq. into Eq. , we obtain

[\frac{\partial}{\partial x} - \frac{i ω}{c_{eff}} - \frac{M_{x} γ_{x}^{2}}{c} v_{⊥} \cdot \nabla_{⊥} - i k_{c} γ_{x}^{2} \sum_{j = 1}^{m} \frac{a_{j, m} (\hat{μ} + \hat{χ})}{1 + b_{j, m} (\hat{μ} + \hat{χ})}] ϕ = 0.

(17)

Equation is the desired WAPE valid for the arbitrary variations in the sound speed and Mach numbers. This characteristic is attained because the Padé series converges for all values of

\hat{μ}

and

\hat{χ}

⁠. Note that Eq. is formulated without the reference wavenumber k₀. This equation is not explicitly presented in Ref. 15. In the Padé (1,1) approximation, only the first term in the series in Eq. is retained. Numerical implementation of this WAPE is considered in Sec. V B. Keeping more terms in the Padé series enables larger propagation angles θ to be considered.

The Padé (1,1) approximation is a more accurate approximation of the square-root operator in Eq. (9) than the first two terms in the Taylor series. (For example, see Table I in Ref. 23.) Therefore, a WAPE in the Padé (1,1) approximation is valid for larger propagation angles than a NAPE. By studying the phase errors pertinent to WAPEs in the Padé (1,1) approximation, Refs. 7 and 15 conclude that these equations can be used if the propagation angles $θ ≲ 35 °$ ⁠. As with NAPEs, this conclusion depends on an acceptable accuracy and the problem at hand.

Section III shows that the dispersion relation obtained from Eq. (17) coincides with geometrical acoustics within the ranges of applicability of Eq. (17). Therefore, Eq. (17) is termed the phase-preserving WAPE. In a stratified moving medium, Eq. (17) coincides with the exact equation, again, within the ranges of its applicability (See Sec. IV.).

As the Padé series converges for all values of the operators ${\hat{μ}}_{σ}$ and ${\hat{χ}}_{σ}$ in Eq. (3), the WAPE given by Eq. (44) in Ref. 15, which is obtained from Eq. (3) with σ = 1, is also phase preserving. However, as explained in Sec. II, more terms in the Padé series might be needed in Eq. (44) from Ref. 15 than in Eq. (17). Furthermore, Eq. (17) is simpler (has fewer terms) than Eq. (44) from Ref. 15. (See also Sec. V B.) Other WAPEs in moving media, such as those in Refs. 2, 13, and 14, are valid for small ε and M.

III. GEOMETRICAL ACOUSTICS

In this section, the ranges of applicability of EWAPEs, NAPEs, and WAPEs are studied by comparing them with the geometrical acoustics equations.

A. Dispersion relation and eikonal equation

In moving media, the dispersion relation for geometrical acoustics [e.g., Eq. (3.35) in Ref. 2] can be written as

k = \frac{ω}{c} - \frac{k \cdot v}{c} .

(18)

Here, k is the wave vector and

k = | k |

⁠. For sound propagation in the positive direction of the x axis, it is convenient to write

k = (k_{x}, k_{⊥})

⁠, where

k_{⊥} = (k_{y}, k_{z})

⁠. Squaring both sides of Eq. and using notations from Eq. yields

(1 - M_{x}^{2}) k_{x}^{2} + 2 M_{x} k_{x} (k_{c} - \frac{k_{⊥} \cdot v_{⊥}}{c}) + k_{⊥}^{2} - k_{c}^{2} + 2 k_{c} \frac{k_{⊥} \cdot v_{⊥}}{c} - \frac{{(k_{⊥} \cdot v_{⊥})}^{2}}{c^{2}} = 0,

(19)

where

k_{⊥} = | k_{⊥} |

⁠. Multiplying both sides of Eq. by

γ_{x}^{2}

⁠, it can be recast in the equivalent form

{[k_{x} + k_{c} M_{x} γ_{x}^{2} (1 - \frac{k_{⊥} \cdot v_{⊥}}{ω})]}^{2} - k_{c}^{2} γ_{x}^{4} [{(1 - \frac{k_{⊥} \cdot v_{⊥}}{ω})}^{2} - \frac{k_{⊥}^{2}}{k_{c}^{2} γ_{x}^{2}}] = 0.

(20)

Equation is an exact consequence of the dispersion relation [Eq. ] and describes sound propagation in the positive and negative directions of the x axis. Equation is recast as

a^{2} - b^{2} = (a - b) (a + b) = 0

⁠, where a² and b² are the first and second terms in this equation, respectively. For sound propagation in the positive direction of the x axis, the latter equation reduces to

a - b = 0

⁠. Thus,

k_{x} + k_{c} M_{x} γ_{x}^{2} (1 - \frac{k_{⊥} \cdot v_{⊥}}{ω}) - k_{c} γ_{x}^{2} \sqrt{{(1 - \frac{k_{⊥} \cdot v_{⊥}}{ω})}^{2} - \frac{k_{⊥}^{2}}{k_{c}^{2} γ_{x}^{2}}} = 0.

(21)

This is the desired dispersion relation for sound propagation in the positive direction of the x axis and is valid for

θ < 90 °

⁠. In Eq. ,

k_{⊥} \sim θ

for small θ. Therefore, the square root has the form

\sqrt{1 + O (θ) + O (θ^{2})}

and can be effectively approximated using the Padé or Taylor series for arbitrary ε and M.

The wave vector k of a sound wave and its phase

Φ

are related by

k = \nabla Φ

⁠. Substitution of this formula into Eq. yields the eikonal equation for the phase

Φ

⁠. For sound propagation in the direction of the x axis, it is convenient to write

(k_{x}, k_{⊥}) = (\frac{\partial Φ}{\partial x}, \nabla_{⊥} Φ) .

(22)

Substituting k_x and

k_{⊥}

from Eq. into Eq. , we can obtain the eikonal equation for the phase

Φ

of a sound wave propagating in the positive direction of the x axis.

B. EWAPE

In this section, starting from Eqs. and , we derive the dispersion relation and compare it with Eq. . The sound pressure p and auxiliary function

ϕ

are expressed in the form used in geometrical acoustics:

p = e^{i k_{0} Θ} p_{0}; p_{0} = p_{1} + \frac{p_{2}}{i k_{0}} + \frac{p_{3}}{{(i k_{0})}^{2}} + \dots,

(23)

ϕ = e^{i k_{0} Θ} ϕ_{0}; ϕ_{0} = ϕ_{1} + \frac{ϕ_{2}}{i k_{0}} + \frac{ϕ_{3}}{{(i k_{0})}^{2}} + \dots .

(24)

Here, Θ is the eikonal and k₀ is the same reference wavenumber as in Eq. . In geometrical acoustics,

1 / k_{0}

plays the role of a small parameter and facilitates formulations of the dispersion relation and transport equation.

Equations and are substituted into Eq. . Keeping terms of order

k_{0}^{(0)}

and omitting higher-order terms, such as

k_{0}^{(- 1)}

⁠, yields

p_{1} = \sqrt{\frac{ϱ}{ϱ_{0}}} (1 - \frac{k_{0} v \cdot \nabla Θ}{ω}) ϕ_{1} .

(25)

The geometrical acoustics equations can be written in forms that do not include the reference wavenumber k₀. To this end, recall that the wave vector of a sound wave can be expressed in terms of the eikonal

k = k_{0} \nabla Θ, (k_{x}, k_{⊥}) = k_{0} (\frac{\partial Θ}{\partial x}, \nabla_{⊥} Θ) .

(26)

Using the first of these relationships, Eq. can be written without k₀ such that

p_{1} = \sqrt{\frac{ϱ}{ϱ_{0}}} (1 - \frac{k \cdot v}{ω}) ϕ_{1} .

(27)

Next, Eq. is substituted into the EWAPE [Eq. ]. To derive the dispersion relation (or the eikonal equation), we keep terms of order k₀ and omit higher-order terms such as

k_{0}^{(0)}, k_{0}^{(- 1)}

⁠, etc. The first term in Eq. is transformed as follows:

\frac{\partial}{\partial x} (e^{i k_{0} Θ} ϕ_{0}) = e^{i k_{0} Θ} (\frac{\partial}{\partial x} + i k_{0} \frac{\partial Θ}{\partial x}) ϕ_{0} = e^{i k_{0} Θ} (\frac{\partial}{\partial x} + i k_{x}) ϕ_{0} = i k_{x} e^{i k_{0} Θ} ϕ_{0} .

(28)

Here, following Eq. ,

k_{0} \partial Θ / \partial x

is replaced with k_x. The last equality is obtained by taking into account that

k_{x} ϕ_{0}

is of order k₀ while

\partial ϕ_{0} / \partial x

is of order

k_{0}^{(0)}

⁠. Thus, when deriving the dispersion relation, the operator

\partial / \partial x

in Eq. is replaced by ik_x, i.e.,

\partial / \partial x \to i k_{x}

⁠. The operator,

{\hat{m}}_{⊥}

⁠, in this equation can be assessed similarly:

{\hat{m}}_{⊥} (e^{i k_{0} Θ} ϕ_{0}) = \frac{i}{ω} e^{i k_{0} Θ} v_{⊥} \cdot (\nabla_{⊥} + i k_{0} \nabla_{⊥} Θ) ϕ_{0} = \frac{i}{ω} e^{i k_{0} Θ} v_{⊥} \cdot (\nabla_{⊥} + i k_{⊥}) ϕ_{0} = m_{⊥} e^{i k_{0} Θ} ϕ_{0}, m_{⊥} = - \frac{k_{⊥} \cdot v_{⊥}}{ω} .

(29)

Here, we used Eq. for

k_{⊥}

and took into account that

k_{⊥} ϕ_{0}

is of order k₀ while

\nabla_{⊥} ϕ_{0}

is of order

k_{0}^{(0)}

⁠. Thus,

{\hat{m}}_{⊥} \to m_{⊥}

⁠.

Equations and exemplify the rule that when one or more differential operators are applied to the ansatz

e^{i k_{0} Θ} ϕ_{0} F

(where F is a function that may appear in transformations), they should act only on the exponential function

e^{i k_{0} Θ}

yielding terms proportional to k₀, whereas the derivatives of F and

ϕ_{0}

can be omitted to order k₀. As a result of this rule,

{\hat{m}}_{⊥}^{n} \to m_{⊥}^{n}

while the operators

{\hat{μ}}_{σ}

and

{\hat{χ}}_{σ}

appearing in Eq. transform to

{\hat{μ}}_{σ} (e^{i k_{0} Θ} ϕ_{0}) = μ_{σ} e^{i k_{0} Θ} ϕ_{0}, μ_{σ} = - \frac{k_{⊥}^{2}}{σ^{2} k_{0}^{2}};

(30)

{\hat{χ}}_{σ} (e^{i k_{0} Θ} ϕ_{0}) = χ_{σ} e^{i k_{0} Θ} ϕ_{0}, χ_{σ} = \frac{k_{c}^{2} γ_{x}^{2}}{k_{0}^{2} σ^{2}} {(1 + m_{⊥})}^{2} - 1.

(31)

To complete the derivation, the square-root operator in Eq. is expressed as the Taylor series and assessed as follows:

\begin{matrix} \sqrt{1 + {\hat{μ}}_{σ} + {\hat{χ}}_{σ}} (e^{i k_{0} Θ} ϕ_{0}) = [1 + \frac{1}{2} ({\hat{μ}}_{σ} + {\hat{χ}}_{σ}) - \frac{1}{8} {({\hat{μ}}_{σ} + {\hat{χ}}_{σ})}^{2} + \dots] (e^{i k_{0} Θ} ϕ_{0}) \\ = [1 + \frac{1}{2} (μ_{σ} + χ_{σ}) - \frac{1}{8} {(μ_{σ} + χ_{σ})}^{2} + \dots] e^{i k_{0} Θ} ϕ_{0} = \sqrt{1 + μ_{σ} + χ_{σ}} e^{i k_{0} Θ} ϕ_{0} . \end{matrix}

(32)

Here, we took into account that the operators

{\hat{μ}}_{σ}

and

{\hat{χ}}_{σ}

commute to order k₀.

Using all these transformations of the operators in Eq. , the dispersion relation following from this equation takes the form

k_{x} + k_{c} M_{x} γ_{x}^{2} (1 + m_{⊥}) - σ k_{0} γ_{x} \sqrt{\frac{k_{c}^{2} γ_{x}^{2}}{k_{0}^{2} σ^{2}} {(1 + m_{⊥})}^{2} + μ_{σ}} = 0.

(33)

Taking into account the definitions of

m_{⊥}

and

μ_{σ}

⁠, it can be shown that this dispersion relation does not depend on σ and coincides with Eq. . Thus, the EWAPE given by Eq. exactly describes the phase of a sound wave within the ranges of its applicability, i.e.,

θ < 90 °

⁠, for any σ.

The EWAPE given by Eq. (3) is derived from the convective Helmholtz equation [Eq. (2)]. The dispersion relation can also be formulated starting with Eq. (2). To this end, Eq. (24) is substituted into the convective Helmholtz equation. Similar to the analysis above, it can be shown that the operator $\nabla$ transforms to $i k$ and the operator Δ transforms to $- k^{2}$ ⁠. With these transformations, the dispersion relation following from Eq. (2) exactly coincides with that used in geometrical acoustics [Eq. (18)].

C. NAPE

This section considers dispersion relations pertinent to NAPEs. To this end, Eq. is substituted into Eq. . In the resulting equation, the operators

v_{⊥} \cdot \nabla_{⊥}

and

Δ_{⊥}

are transformed similarly to Eqs. and , respectively. As a result, we obtain the following dispersion relation:

k_{x} + \frac{k_{⊥}^{2}}{2 k_{c}} - \frac{ω}{c_{eff}} + \frac{k_{⊥} \cdot v_{⊥}}{c_{eff}} = 0.

(34)

The same result can be derived from the dispersion relation in geometrical acoustics [Eq. ]. To this end, the square root in Eq. should be approximated similarly to Eq. . Thus, within the ranges of its applicability, the NAPE [Eq. ] correctly describes the phase of a sound wave for arbitrary variations in the sound speed and Mach numbers.

Other NAPEs in moving media correctly describe the phase of a sound wave only approximately when M and/or ε are small. For example, starting with Eq. (14), the dispersion relation can be derived. It can be shown that it coincides with the dispersion relation given by Eq. (34) only for small ε and M and the NAPE given by Eq. (46) from Ref. 15 preserves the phase of a sound wave only for small ε.

D. WAPE

Next, we consider the dispersion relation obtained from the WAPE given by Eq. , where the auxiliary function

ϕ

is specified by Eq. . The operators

\hat{μ}

and

\hat{χ}

are estimated with Eqs. and , respectively, in which σ is given by Eq. such that

\hat{μ} (e^{i k_{0} Θ} ϕ_{0}) = μ e^{i k_{0} Θ} ϕ_{0}, μ = - \frac{k_{⊥}^{2}}{k_{c}^{2} γ_{x}^{2}};

(35)

\hat{χ} (e^{i k_{0} Θ} ϕ_{0}) = χ e^{i k_{0} Θ} ϕ_{0}, χ = 2 m_{⊥} + m_{⊥}^{2} .

(36)

Using these results, Eq. becomes

k_{x} - \frac{ω}{c_{eff}} - \frac{M_{x} γ_{x}^{2}}{c} k_{⊥} \cdot v_{⊥} - k_{c} γ_{x}^{2} \sum_{j = 1}^{m} \frac{a_{j, m} (μ + χ)}{1 + b_{j, m} (μ + χ)} = 0.

(37)

The same equation can be obtained from the dispersion relation in geometrical acoustics [Eq. ]. To this end, the square root in Eq. is expressed as the Padé series [Eq. ]. Therefore, in the Padé approximation, the phase-preserving WAPE [Eq. ] correctly describes the phase of a sound wave for arbitrary variations in the sound speed and Mach numbers, i.e., arbitrary ε and M.

The WAPE obtained from Eq. (3) in the Padé approximation also preserves the phase of a sound wave within the ranges of its applicability but might require more terms in the corresponding Padé series. (See Sec. II.) Other WAPEs appearing in the literature^2,13,14 preserve the phase of a sound wave only for small ε and M.

IV. STRATIFIED MOVING MEDIUM

In this section, the ranges of applicability of EWAPEs, NAPEs, and WAPEs are studied by comparing them to an exact result for sound propagation in a stratified moving medium.

A. Exact equation

As is commonly done in atmospheric acoustics, here, we assume a stratified moving medium. Thus, the ambient quantities depend only on the vertical coordinate z, i.e.,

c = c (z), ϱ = ϱ (z)

⁠, v_z = 0, and

{\overset{ˇ}{v}}_{⊥} (z) = (v_{x} (z), v_{y} (z))

⁠. Equations (2.63) and (2.64) in Ref. 2 exactly describe sound propagation in a stratified moving medium. For a monochromatic sound wave, by omitting the time dependence

e^{- i ω t}

⁠, these equations can be recast as

p (\overset{ˇ}{r}, z) = \sqrt{\frac{ϱ}{ϱ_{0}}} \int_{- \infty}^{\infty} (1 - \frac{κ \cdot {\overset{ˇ}{v}}_{⊥}}{ω}) \hat{ϕ} (κ, z) exp (i κ \cdot \overset{ˇ}{r}) d^{2} κ,

(38)

[\frac{\partial^{2}}{\partial z^{2}} + \frac{{(ω - κ \cdot {\overset{ˇ}{v}}_{⊥})}^{2}}{c^{2}} - κ^{2}] \hat{ϕ} (κ, z) = 0.

(39)

Here,

\overset{ˇ}{r} = (x, y), κ = (κ_{x}, κ_{y})

⁠, and the function

\hat{ϕ} (κ, z)

is proportional to the two-dimensional (2D) spatial spectrum of the sound pressure. In this section, the check marks above

\overset{ˇ}{r}

and

{\overset{ˇ}{v}}_{⊥}

are used to distinguish these quantities from those in Secs. II and III. In Eq. , the derivatives of c, ϱ, and

{\overset{ˇ}{v}}_{⊥}

are omitted, as consistent with the high-frequency approximation considered in this article.

Equation can be recast in the following equivalent form:

[(1 - M_{x}^{2}) κ_{x}^{2} + 2 k_{c} M_{x} (1 - \frac{κ_{y} v_{y}}{ω}) κ_{x} - \frac{\partial^{2}}{\partial z^{2}} + κ_{y}^{2} - k_{c}^{2} {(1 - \frac{κ_{y} v_{y}}{ω})}^{2}] \hat{ϕ} = 0.

(40)

Multiplying both sides of this equation by

γ_{x}^{2}

⁠, it can be written as

{{[κ_{x} + k_{c} M_{x} γ_{x}^{2} (1 - \frac{κ_{y} v_{y}}{ω})]}^{2} - k_{c}^{2} γ_{x}^{4} [{(1 - \frac{κ_{y} v_{y}}{ω})}^{2} + \frac{1}{k_{c}^{2} γ_{x}^{2}} \frac{\partial^{2}}{\partial z^{2}} - \frac{κ_{y}^{2}}{k_{c}^{2} γ_{x}^{2}}]} \hat{ϕ} = 0.

(41)

Because we omit the derivatives of c and

{\overset{ˇ}{v}}_{⊥}

⁠, Eq. can be recast as

a^{2} - b^{2} = (a - b) (a + b) = 0

⁠. For sound propagation in the positive direction of the x axis, the latter equation reduces to

a - b = 0

⁠. Thus,

[κ_{x} + k_{c} M_{x} γ_{x}^{2} (1 - \frac{κ_{y} v_{y}}{ω}) - k_{c} γ_{x}^{2} \sqrt{{(1 - \frac{κ_{y} v_{y}}{ω})}^{2} + \frac{1}{k_{c}^{2} γ_{x}^{2}} \frac{\partial^{2}}{\partial z^{2}} - \frac{κ_{y}^{2}}{k_{c}^{2} γ_{x}^{2}}}] \hat{ϕ} = 0.

(42)

Equation exactly describes sound propagation in the positive direction of the x axis (i.e.,

θ < 90 °

⁠) provided that the derivatives of c and

{\overset{ˇ}{v}}_{⊥}

can be omitted.

B. EWAPE

Next, we consider the EWAPE [Eq. ] in a stratified moving medium. The auxiliary function

ϕ

is expressed in the form

ϕ (\overset{ˇ}{r}, z) = \int \hat{ϕ} (κ, z) exp (i κ \cdot \overset{ˇ}{r}) d^{2} κ .

(43)

Substitution of Eq. into Eq. yields Eq. . Thus, the relationship between p and

\hat{ϕ}

is the same as that in Ref. 2 and the formulations pertinent to the EWAPE given by Eq. . When Eq. is substituted into Eq. , the operators are transformed as follows:

\frac{\partial}{\partial x} \to i κ_{x},

(44)

{\hat{m}}_{⊥} = \frac{i}{ω} {\overset{ˇ}{v}}_{⊥} \cdot \nabla_{⊥} = \frac{i v_{y}}{ω} \frac{\partial}{\partial y} \to m_{⊥, s} = - \frac{κ_{y} v_{y}}{ω},

(45)

{\hat{μ}}_{σ} = \frac{1}{σ^{2} k_{0}^{2}} (\frac{\partial^{2}}{\partial y^{2}} + \frac{\partial^{2}}{\partial z^{2}}) \to {\hat{μ}}_{σ, s} = \frac{1}{σ^{2} k_{0}^{2}} (\frac{\partial^{2}}{\partial z^{2}} - κ_{y}^{2}),

(46)

{\hat{χ}}_{σ} \to χ_{σ, s} = \frac{k_{c}^{2} γ_{x}^{2}}{σ^{2} k_{0}^{2}} {(1 + m_{⊥, s})}^{2} - 1.

(47)

Here, the subscript s indicates that the corresponding quantity is pertinent to a stratified medium. When deriving these relationships, we also took into account that the operators

{\hat{m}}_{⊥}, {\hat{μ}}_{σ}

⁠, and

{\hat{χ}}_{σ}

commute if the derivatives of c and

{\overset{ˇ}{v}}_{⊥}

are omitted. Because of this property, in Eq. ,

\sqrt{1 + {\hat{μ}}_{σ} + {\hat{χ}}_{σ}} \to \sqrt{1 + {\hat{μ}}_{σ, s} + χ_{σ, s}} .

(48)

With these results, Eq. for

ϕ

transforms to the following equation for

\hat{ϕ}

⁠:

[κ_{x} + k_{c} M_{x} γ_{x}^{2} (1 + m_{⊥, s}) - σ k_{0} γ_{x} \sqrt{\frac{k_{c}^{2} γ_{x}^{2}}{σ^{2} k_{0}^{2}} {(1 + m_{⊥, s})}^{2} + \frac{1}{σ^{2} k_{0}^{2}} \frac{\partial^{2}}{\partial z^{2}} - \frac{κ_{y}^{2}}{σ^{2} k_{0}^{2}}}] \hat{ϕ} = 0.

(49)

This is the EWAPE for propagation in a stratified moving medium, which is derived from Eq. . Substituting with

m_{⊥, s}

⁠, it can be shown that Eq. does not depend on σ and coincides with Eq. . Thus, the EWAPE given by Eq. exactly describes propagation in a stratified moving medium within the ranges of its applicability (for

θ < 90 °

and when the derivatives of c and

{\overset{ˇ}{v}}_{⊥}

can be omitted) for any σ.

C. NAPE

Formulations for the NAPE proceed similarly to those for the EWAPE. To this end, Eq. is substituted into Eq. . The operators

\partial / \partial x

and

{\overset{ˇ}{v}}_{⊥} \cdot \nabla_{⊥}

are assessed with Eqs. and . The result is

[κ_{x} - \frac{ω}{c_{eff}} - \frac{1}{2 k_{c}} (\frac{\partial^{2}}{\partial z^{2}} - κ_{y}^{2}) + \frac{κ_{y} v_{y}}{c_{eff}}] \hat{ϕ} (κ, z) = 0.

(50)

This is the NAPE for a stratified moving medium, which follows from Eq. . The same equation can be obtained from Eq. by approximating the square root with Eq. .

Thus, within the ranges of its applicability, the phase-preserving NAPE [Eq. (13)] correctly describes sound propagation in a stratified moving medium with arbitrary variations in the sound speed ε and Mach number M. Other NAPEs [such as Eq. (14)] are valid only if M and/or ε are small.

D. WAPE

Next, consider the WAPE [Eq. ] in a stratified moving medium, where

ϕ

is specified by Eq. . To derive the equation for

\hat{ϕ}

⁠, the operators

\partial / \partial x

and

{\overset{ˇ}{v}}_{⊥} \cdot \nabla_{⊥}

are replaced with Eqs. and . With σ given by Eq. , Eqs. and for the operators

{\hat{μ}}_{σ}

and

{\hat{χ}}_{σ}

⁠, respectively, take the form

\hat{μ} \to {\hat{μ}}_{s} = \frac{1}{k_{c}^{2} γ_{x}^{2}} (\frac{\partial^{2}}{\partial z^{2}} - κ_{y}^{2}), \hat{χ} \to χ_{s} = 2 m_{⊥, s} + m_{⊥, s}^{2} .

(51)

Here, the subscript σ is omitted because the factor σ is specified with Eq. . Also, similar to Eq. , the following transformation holds:

\frac{a_{j, m} (\hat{μ} + \hat{χ})}{1 + b_{j, m} (\hat{μ} + \hat{χ})} \to \frac{a_{j, m} ({\hat{μ}}_{s} + χ_{s})}{1 + b_{j, m} ({\hat{μ}}_{s} + χ_{s})} .

(52)

With these results, the WAPE for

\hat{ϕ}

takes the form

[κ_{x} - \frac{ω}{c_{eff}} - M_{x} γ_{x}^{2} \frac{κ_{y} v_{y}}{c} - k_{c} γ_{x}^{2} \sum_{j = 1}^{m} \frac{a_{j, m} ({\hat{μ}}_{s} + χ_{s})}{1 + b_{j, m} ({\hat{μ}}_{s} + χ_{s})}] \hat{ϕ} (κ, z) = 0.

(53)

This is the WAPE in a stratified moving medium, derived from Eq. . The same equation can be obtained from Eq. by approximating the square root with the Padé series [Eq. ]. Thus, within the ranges of its applicability, the phase-preserving WAPE [Eq. ] coincides with the exact equation for sound propagation in a stratified moving medium for arbitrary variations in the sound speed ε and Mach numbers M.

It can be shown that the WAPE given by Eq. (44) in Ref. 15 also coincides with the exact equation for sound propagation in a stratified medium within the ranges of its applicability. However, this equation might require more terms in the Padé series than in Eq. (53), as discussed in Sec. II. The WAPEs in Refs. 2, 13, and 14 coincide with the exact equation in a stratified moving medium only approximately for small ε and M.

V. NUMERICAL IMPLEMENTATION OF THE PARABOLIC EQUATIONS

In this section, we consider numerical implementation of the phase-preserving NAPE [Eq. (13)] and WAPE [Eq. (17)] for 2D sound propagation in the vertical (x,z) plane. As is typical for sound propagation modeling, we also assume that v_z = 0.

A. NAPE

First, consider the phase-preserving NAPE. For 2D sound propagation, the operators in Eq. simplify such that

v_{⊥} \cdot \nabla_{⊥} = 0

and

Δ_{⊥} = \partial^{2} / \partial z^{2}

⁠. In this equation, the auxiliary function

ϕ

is expressed in the form

ϕ (x, z) = \hat{ϕ} (x, z) e^{i k_{0} x} .

(54)

Here, k₀ is the reference wavenumber and

\hat{ϕ}

is the complex amplitude of the auxiliary function. With these results, the phase-preserving NAPE [Eq. ] takes the form

[\frac{\partial}{\partial x} - \frac{i}{2 k_{0} n} \frac{\partial^{2}}{\partial z^{2}} - i k_{0} (n_{eff} - 1)] \hat{ϕ} = 0,

(55)

where n is the refractive index in a motionless medium and

n_{eff}

is the effective refractive index,

n = \frac{c_{0}}{c}, n_{eff} = \frac{c_{0}}{c_{eff}} = \frac{c_{0}}{c + v_{x}} .

(56)

Note that

n_{eff}

appears in Eq. as a convenient notation rather than the result of the effective sound speed approximation.

For further analysis, it is convenient to express Eq. in the equivalent form

(h_{1, 0} + \frac{h_{1, 2}}{k_{0}^{2}} \frac{\partial^{2}}{\partial z^{2}}) \frac{\partial \hat{ϕ}}{\partial x} = i k_{0} (h_{2, 0} + \frac{h_{2, 2}}{k_{0}^{2}} \frac{\partial^{2}}{\partial z^{2}}) \hat{ϕ},

(57)

where the coefficients

h_{n, m}

are given by

h_{1, 0} = 1, h_{1, 2} = 0, h_{2, 0} = n_{eff} - 1, h_{2, 2} = \frac{1}{2 n} .

(58)

The initial condition for the function

\hat{ϕ}

at x = 0 is chosen as that adopted in the literature for a point source and NAPE,^1,6

\hat{ϕ} (x = 0, z) = - \frac{i ϕ_{0}}{2 \sqrt{2 π}} exp (- \frac{k_{0}^{2} {(z - z_{s})}^{2}}{2}) .

(59)

Here,

ϕ_{0}

characterizes the amplitude of the point source and z_s is the source height above the ground. In numerical examples to follow,

ϕ_{0} = 1

Pa. Numerical implementation of Eq. is outlined in Sec. VI A of Ref. 15. It involves only minimal modifications to existing Crank-Nicholson NAPE algorithms and is just as computationally efficient.

After the auxiliary function

\hat{ϕ}

is found as a solution of Eq. , the sound pressure is calculated as

p (x, z) = \sqrt{\frac{ϱ}{ϱ_{0}}} e^{i k_{0} x} (1 - M_{x} + \frac{i M_{x}}{k_{0}} \frac{\partial}{\partial x}) \hat{ϕ} (x, z) .

(60)

This result follows from Eqs. and . Numerical implementation of Eq. is also considered in Sec. VI A of Ref. 15.

In the literature, the sound pressure in a moving medium, such as the atmosphere, is often calculated with the NAPE based on the effective sound speed approximation or the NAPE_eff.¹ This equation is given by Eq. with the last term omitted. In the NAPE_eff, the sound pressure is expressed in the form

p (x, z) = \hat{p} (x, z) e^{i k_{0} x},

(61)

where

\hat{p}

is the complex amplitude of the sound pressure. For 2D sound propagation calculations using Eqs. and ,

\hat{p}

also satisfies Eq. but with different coefficients

h_{n, m}

⁠, namely,

h_{1, 0} = 1, h_{1, 2} = 0, h_{2, 0} = \frac{n_{eff}^{2} - 1}{2}, h_{2, 2} = \frac{1}{2} .

(62)

The initial condition for the NAPE_eff is formulated for

\hat{p} (0, z)

and is given by the right-hand side of Eq. .

B. WAPE

Next, consider the phase-preserving WAPE [Eq. ]. For 2D sound propagation in the (x,z)-plane, the operators

\hat{χ}

and

\hat{μ}

in this equation simplify to

\hat{χ} = 0, \hat{μ} = \frac{1}{k_{0}^{2} n^{2} γ_{x}^{2}} \frac{\partial^{2}}{\partial z^{2}} .

(63)

In Eq. , ϕ is replaced with Eq. . With these transformations and in the Padé (1,1) approximation, the phase-preserving WAPE takes the form

[\frac{\partial}{\partial x} - i k_{0} (n_{eff} - 1) - i k_{0} n γ_{x}^{2} \frac{a_{1, 1} \hat{μ}}{1 + b_{1, 1} \hat{μ}}] \hat{ϕ} = 0.

(64)

It follows from Eq. that in the Padé (1,1) approximation,

a_{1, 1} = 1 / 2

and

b_{1, 1} = 1 / 4

⁠. Multiplying both sides of Eq. by

1 + b_{1, 1} \hat{μ}

⁠, this equation can be recast as Eq. , in which the coefficients

h_{n, m}

are given by

h_{1, 0} = 1, h_{1, 2} = \frac{1}{4 n^{2} γ_{x}^{2}}, h_{2, 0} = n_{eff} - 1, h_{2, 2} = \frac{1}{2 n} + \frac{n_{eff} - 1}{4 n^{2} γ_{x}^{2}} .

(65)

Although written in a somewhat different form, Eq. for the coefficients

h_{n, m}

is consistent with Eq. (B5) in Ref. 15. The initial condition for the WAPE in the Padé (1,1) approximation is chosen as that adopted in the literature [see Eq. (G.75) in Ref. 1]:

\hat{ϕ} (x = 0, z) = - \frac{i ϕ_{0}}{2 \sqrt{2 π}} [1.3717 - 0.3701 k_{0}^{2} {(z - z_{s})}^{2}] exp (- \frac{k_{0}^{2} {(z - z_{s})}^{2}}{3}) .

(66)

In the literature, the sound pressure p in a moving medium is often calculated with the WAPE based on the effective sound speed approximation or the WAPE_eff, e.g., see Ref. 1 and Eq. (56) in Ref. 15. In the Padé (1,1) approximation and substituting with Eq. , the WAPE_eff can also be written as Eq. with the following coefficients:

h_{1, 0} = 1 + \frac{n_{eff}^{2} - 1}{4}, h_{1, 2} = \frac{1}{4}, h_{2, 0} = \frac{n_{eff}^{2} - 1}{2}, h_{2, 2} = \frac{1}{2} .

(67)

The initial condition for the WAPE_eff coincides with the right-hand side of Eq. .

C. Intercomparison of the parabolic equations

It follows from Secs. V A and V B that the phase-preserving NAPE and WAPE and the corresponding equations based on the effective sound speed approximation, NAPE_eff and WAPE_eff, reduce to the same equation [Eq. (57)] but with different coefficients $h_{n, m}$ ⁠. These coefficients are summarized in Table I. The rows “NAPE” and “WAPE” correspond to the phase-preserving NAPE and WAPE [Eqs. (13) and (17)], and the rows “NAPE_eff” and “WAPE_eff” correspond to the equations based on the effective sound speed approximation. The row “WAPE_M” corresponds to Eq. (78) in Ref. 15, which is the 2D WAPE in the Padé (1,1) approximation and can be obtained from Eq. (3) in this article by setting σ = 1 and written in the form given by Eq. (57). As mentioned in the Introduction, Eq. (78) from Ref. 15 has already been used in the literature^16–19 and is valid for arbitrary Mach numbers but small variations in the sound speed. In the WAPE_M row, $τ = (n - M_{x}) M_{x} γ_{x}^{2}$ ⁠.

Interestingly, although the parabolic equations in Table I describe similar sound propagation phenomena, the coefficients $h_{n, m}$ in these equations differ significantly. Also, these coefficients are more complex for the WAPE_M than for the phase-preserving WAPE. This result elucidates the statements above that the EWAPE, NAPE, and WAPE, which are obtained from Eq. (3) by setting σ equal to Eq. (8), are simpler than those pertinent to σ = 1.

Another difference between these parabolic equations is that the phase-preserving NAPE and WAPE and the WAPE_M are formulated for the complex amplitude of the auxiliary function $\hat{ϕ}$ ⁠, while the NAPE_eff and WAPE_eff are formulated for the complex amplitude of the sound pressure $\hat{p}$ ⁠.

VI. NUMERICAL RESULTS

In this section, numerical results obtained with the phase-preserving NAPE and WAPE and the NAPE_eff and WAPE_eff are compared to analytical and reference solutions. We will also briefly consider results obtained with the WAPE_M.

In all numerical examples to follow, we consider a point source with the frequency $f = ω / (2 π) = 200$ Hz located above a perfectly reflecting, flat surface. The reference sound speed is $c_{0} = 340$ m/s. The computational domain is 500 m long and 125 m high. At the lower boundary, the vertical derivative of $ϕ$ for NAPE and WAPE and the vertical derivative of p for NAPE_eff and WAPE_eff are zero due to a perfectly reflecting surface. For the upper boundary, a perfectly matched layer (PML) is used.²⁶ The thickness of the PML is 30 grid points. It is implemented as described in Ref. 19. In all numerical examples, for x greater than about 100 m but smaller than 500 m, the propagation angles $θ ≲ 35 °$ ⁠.

A. Uniformly moving medium

This section considers the point source located 50 m above the surface in a uniformly moving medium with Mach number M = 0.5 and constant sound speed c = c₀. The relatively large Mach number is chosen to elucidate the difference between the results obtained with the phase-preserving NAPE and WAPE and the NAPE_eff and WAPE_eff. Note that for the considered case c = c₀, the WAPE_M coincides with the phase-preserving WAPE, i.e., the coefficients h_nm in Table I for these parabolic equations are the same. Using these parabolic equations, the sound pressure level (SPL) and relative level are calculated. The latter is defined as follows

Δ L (x, z) = 20 {log}_{10} (\frac{| p (x, z) |}{| p_{0} (x, z) |}),

(68)

where p₀ is the sound pressure produced by the point source in free space. In these calculations, the steps along the x and z axes are 0.05 m, which corresponds to

λ / 34

⁠. Here, λ is the sound wavelength. These relatively small steps are needed to correctly predict the sound pressure at interference minima. In practical applications, the grid steps in the WAPE and NAPE can be set to those typically used in the literature, e.g.,

λ / 10

⁠.¹ Therefore, the computational cost of these equations is about the same as for the NAPE_eff and WAPE_eff.

The SPL calculated with the parabolic equations is depicted in Fig. 1, which also presents the SPL corresponding to the exact analytical solution of the problem considered, presented in Sec. VI B of Ref. 15. The SPL has several minima and maxima due to interference between the direct wave from the source to receiver and that reflected from the surface. For $x ≳ 100$ m, the SPL calculated with the phase-preserving WAPE is indistinguishable from the exact analytical solution. The phase-preserving NAPE results are also close to the analytical solution but at larger propagation ranges. The results obtained with the WAPE_eff and NAPE_eff noticeably deviate from the analytical solution.

FIG. 1.

View large Download slide

(Color online) SPL in a uniformly moving medium with M = 0.5 produced by a point source located 50 m above a perfectly reflecting, flat surface. The sound frequency is 200 Hz. The subplots correspond to the exact analytical solution and numerical results obtained with the parabolic equations as labeled.

These conclusions are reinforced by Fig. 2, in which the relative level $Δ L$ obtained with the analytical solution and parabolic equations, is plotted versus the propagation range x for four heights z above the surface. For $x ≲ 100$ m, propagation angles θ can be relatively large such that the ranges of applicability of the NAPEs and WAPEs might not be fulfilled. Therefore, the predictions based on the parabolic equations can significantly deviate from the analytical solution. For $x ≳ 100$ m, the results obtained with the phase-preserving WAPE are very close to the analytical solution. The results obtained with the phase-preserving NAPE slightly deviate from the analytical solution in the range $100 m ≲ x ≲ 400$ m because the NAPE is valid for smaller θ than the WAPE. However, for $x ≳ 400$ m (where θ becomes small), the NAPE is close to the analytical solution. On the other hand, Fig. 2 clearly shows that the results obtained with the WAPE_eff and NAPE_eff significantly deviate from the analytical solution. For the problem considered (but with different source height and frequency), inaccurate predictions of WAPE_eff are also reported in Ref. 15.

FIG. 2.

View large Download slide

(Color online) Relative level $Δ L$ versus the propagation range x for four heights z above the surface. The geometry of the problem is the same as that for the SPL in Fig. 1. Different curves correspond to the analytical solution and parabolic equations.

The difference between the results obtained with the parabolic equations and analytical solution can be quantified by the normalized errors,

ϵ (x, z) = \frac{| p_{PE} (x, z) | - | p_{an} (x, z) |}{| p_{0} (x, z) |},

(69)

where

p_{PE}

is the sound pressure calculated with one of the parabolic equations from Table I, and

p_{an}

is the sound pressure obtained with the exact analytical solution. (Note that the errors ϵ are normalized by p₀ rather than by

p_{an}

because

p_{an}

can be very small at the interference minima; see Fig. 2.) The normalized errors ϵ versus the propagation range x are plotted in Fig. 3 for four heights z above the ground. These errors correspond to the results shown in Fig. 2. In Fig. 3, the x axis starts with

x_{min} = 100

m because at shorter ranges, the parabolic equations might not applicable. It follows from Fig. 3 that the normalized errors ϵ pertinent to the phase-preserving WAPE are almost zero. The normalized errors ϵ for the phase-preserving NAPE are small and tend to zero as the propagation range x increases. The normalized errors pertinent to the WAPE_eff and NAPE_eff are large and, generally, do not decrease with range.

FIG. 3.

View large Download slide

(Color online) Normalized errors ϵ between the analytical solution and results obtained with the parabolic equations from Fig. 2.

Next, we consider the cumulative normalized errors, which are obtained by integrating the normalized errors ϵ along the x axis:

\bar{ϵ} (z) = \frac{1}{x_{max} - x_{min}} \int_{x_{min}}^{x_{max}} | ϵ (x, z) | d x,

(70)

where

x_{max} = 500

m. The cumulative errors pertinent to the normalized errors in Fig. 3 are shown in Table II. It follows from the table that the phase-preserving WAPE has the smallest cumulative errors, followed by those for the phase-preserving NAPE and the NAPE_eff and WAPE_eff.

TABLE II.

Cumulative errors $\bar{ϵ}$ pertinent to the normalized errors ϵ shown in Fig. 3.

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.03	0.03	0.09	0.15
WAPE	0.00	0.01	0.01	0.02
NAPE_eff	0.44	0.32	0.45	0.43
WAPE_eff	0.46	0.51	0.59	0.76

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.03	0.03	0.09	0.15
WAPE	0.00	0.01	0.01	0.02
NAPE_eff	0.44	0.32	0.45	0.43
WAPE_eff	0.46	0.51	0.59	0.76

B. Linear sound speed profile

This section considers sound propagation in a motionless medium. The source is located 10 m above the surface. The sound speed linearly increases with the height z such that

c (z) = c_{0} + α z,

(71)

where α = 1 s⁻¹. As in Sec. VI A, the relatively large value of the coefficient α is chosen to elucidate the difference between the results obtained with the phase-preserving WAPE and NAPE and the WAPE_eff and NAPE_eff. It follows from Table I that in a motionless medium, the WAPE_M coincides with the WAPE_eff.

The SPL calculated with these parabolic equations is depicted in Fig. 4, which also shows the SPL calculated with the finite-difference time-domain (FDTD) solution of linearized equations of fluid dynamics. (For details, see the Appendix.) For the linear sound speed profile, the FDTD approach can be considered as a reference solution.

FIG. 4.

View large Download slide

(Color online) SPL for a linear profile of the sound speed produced by a point source located 10 m above a perfectly reflecting, flat surface. The sound frequency is 200 Hz. The subplots correspond to the results obtained with the FDTD solution and parabolic equations as labeled.

It can be observed from Fig. 4 that the phase-preserving WAPE provides the best agreement with the FDTD solution, followed by the WAPE_eff, NAPE, and NAPE_eff. This conclusion is further illustrated in Fig. 5, where the relative level $Δ L$ calculated with the FDTD approach and parabolic equations is plotted versus the propagation range x for four heights z above the ground. Figure 5 clearly shows that for $x ≳ 100$ m, the results obtained with the phase-preserving WAPE are in a very good agreement with those obtained with the FDTD solution.

FIG. 5.

View large Download slide

(Color online) Relative level $Δ L$ versus the propagation range x for a linear profile of the sound speed and four heights z above the surface. The geometry of the problem is the same as that for the SPL in Fig. 4 Different curves correspond to the results obtained with the FDTD solution and parabolic equations.

The normalized errors ϵ pertinent to the results displayed in Fig. 5 are depicted in Fig. 6. It follows from Fig. 6 that the normalized errors for the phase-preserving WAPE are smaller compared to those for the WAPE_eff. Similarly, the phase-preserving NAPE has smaller normalized errors compared to the NAPE_eff. The WAPE_eff outperforms the NAPE because this section considers sound propagation in a motionless medium for which the effective sound speed approximation is not relevant.

FIG. 6.

View large Download slide

(Color online) Normalized errors ϵ between the FDTD solution and results obtained with the parabolic equations from Fig. 5.

Table III shows the cumulative errors $\bar{ϵ}$ pertinent to the normalized errors ϵ in Fig. 6. It follows from Table III that the phase-preserving WAPE has smaller cumulative errors compared to the WAPE_eff, as does the phase-preserving NAPE compared to the NAPE_eff.

TABLE III.

Cumulative errors $\bar{ϵ}$ pertinent to the normalized errors ϵ depicted in Fig. 6.

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.34	0.51	0.68	1.01
WAPE	0.08	0.09	0.11	0.17
NAPE_eff	0.52	1.09	1.15	1.55
WAPE_eff	0.27	0.41	0.50	0.69

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.34	0.51	0.68	1.01
WAPE	0.08	0.09	0.11	0.17
NAPE_eff	0.52	1.09	1.15	1.55
WAPE_eff	0.27	0.41	0.50	0.69

C. Stratified moving medium

This section considers a stratification for which the sound speed and medium velocity depend on z. Again, the sound speed profile is given by Eq. but with a much smaller coefficient

α = 0.1

s⁻¹. The medium velocity logarithmically increases with height,

v_{x} (z) = v_{0} ln (1 + β z),

(72)

where

v_{0} = 9

m/s and β = 1 m⁻¹. As in Sec. VI B, the source is located at z_s = 10 m above the surface.

Figure 7 shows the relative level $Δ L$ versus the propagation range x for four heights z above the surface. The results correspond to the FDTD solution, the phase-preserving NAPE and WAPE, and the NAPE_eff and WAPE_eff. For the considered case of a stratified moving medium, the FDTD solution is valid for arbitrary Mach numbers (see Sec. IV D in Ref. 27). It follows from Fig. 7 that the results obtained with the phase-preserving WAPE are very close to the FDTD solution and the phase-preserving NAPE outperforms the NAPE_eff and WAPE_eff. The results obtained with the latter two parabolic equations often deviate from the reference solution not only at the interference minima and maxima but also at the ranges between them.

FIG. 7.

View large Download slide

(Color online) Relative level $Δ L$ versus the propagation range x for a stratified moving medium (Sec. V C) and four heights z above a perfectly reflecting, flat surface. The source is located 10 m above the surface and the sound frequency is 200 Hz. Different curves correspond to the results obtained with the FDTD solution and parabolic equations.

Table IV shows the cumulative errors $\bar{ϵ}$ pertinent to the results in Fig. 7. It follows from the table that the phase-preserving WAPE has the smallest errors, followed by the phase-preserving NAPE, WAPE_eff, and NAPE_eff. The WAPE_M was also used to calculate the relative level for the stratification considered in this section. The results obtained are very close to those for the phase-preserving WAPE with the same cumulative errors. This can be explained by the fact that for the stratification considered, the variation in the sound speed is relatively small.

TABLE IV.

Cumulative errors $\bar{ϵ}$ pertinent to the results depicted in Fig. 7.

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.09	0.08	0.13	0.13
WAPE	0.04	0.02	0.02	0.05
NAPE_eff	0.18	0.80	0.56	0.92
WAPE_eff	0.05	0.14	0.17	0.21

Equation	z = 1 m	z = 2 m	z = 5 m	z = 10 m
NAPE	0.09	0.08	0.13	0.13
WAPE	0.04	0.02	0.02	0.05
NAPE_eff	0.18	0.80	0.56	0.92
WAPE_eff	0.05	0.14	0.17	0.21

VII. CONCLUSIONS

This article considered the phase-preserving EWAPEs, WAPEs, and NAPEs for sound propagation in moving media. These equations [Eqs. (9), (17), and (13), respectively] were derived in the high-frequency approximation while omitting the derivatives of the sound speed and medium velocity. Within the ranges of their applicability, the EWAPE, WAPE, and NAPE exactly describe the phase of sound waves and are valid for arbitrary variations in the sound speed and arbitrary (subsonic) Mach numbers. These equations also correctly describe sound propagation in a stratified moving medium within the ranges of their applicability.

Although variations in sound speed and wind Mach number are often relatively small, omitting the corresponding second-order terms can result in significant phase errors for long propagation ranges. Therefore, it is preferable to use the phase-preserving WAPE and NAPE, at least in some practical applications.

This article further showed that the WAPE given by Eq. (44) in Ref. 15 can also correctly describe the phase of a sound wave but might require more terms in the Padé series expansion. Moreover, that equation involves more terms than the phase-preserving WAPE given by Eq. (17) in this article. The WAPEs from Refs. 2, 13, and 14 preserve the phase of sound waves only approximately for small variations in sound speed and Mach numbers. Previous NAPEs from the literature (e.g., Refs. 2 and 15) correctly describe the phase of sound waves only for small variations in the sound speed and/or low Mach numbers.

Numerical simulations showed that the results obtained with the phase-preserving WAPE and NAPE agree better with the analytical and reference solutions than those with the parabolic equations based on the effective sound speed approximation (the WAPE_eff and NAPE_eff), which are often used in the literature. Limitations of this approximation have been reported elsewhere, e.g., Refs. 2, 15, and 28. There are, of course, problems for which the effective sound speed approximation might be suitable. However, because the complexity of numerical implementation of the phase-preserving WAPE and NAPE is about the same as that for the WAPE_eff and NAPE_eff, it makes sense to use the phase-preserving equations.

ACKNOWLEDGMENTS

This research was partially funded by the U.S. Army Engineer Research and Development Center (ERDC) basic research program. Permission to publish was granted by Director, Cold Regions Research and Engineering Laboratory. Any opinions expressed in this article are those of the authors and are not to be construed as official positions of the funding agency or the U.S. Department of the Army unless so designated by other authorized documents. This work was also performed within the framework of the LABEX CeLya (Grant No. ANR-10-LABX-0060) of Université de Lyon within the program Investissements d'Avenir (Grant No. ANR-16-IDEX-0005), operated by the French National Research Agency (ANR).

AUTHOR DECLARATION

Conflict of Interest

The authors have no conflicts to disclose.

DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this study.

Appendix

This appendix briefly outlines the FDTD approach employed in Secs. VI B and VI C. The approach is based on solution of Eqs. (17) and (18) in Ref. 27, which, for 2D sound propagation in a stratified moving medium, become

\begin{matrix} \frac{\partial p}{\partial t} + v_{x} \frac{\partial p}{\partial x} + ϱ c^{2} (\frac{\partial w_{x}}{\partial x} + \frac{\partial w_{z}}{\partial z}) = ϱ c^{2} Q, \\ \frac{\partial w_{x}}{\partial t} + v_{x} \frac{\partial w_{x}}{\partial x} + \frac{1}{ϱ} \frac{\partial p}{\partial x} = 0, \\ \frac{\partial w_{z}}{\partial t} + v_{x} \frac{\partial w_{z}}{\partial x} + \frac{1}{ϱ} \frac{\partial p}{\partial z} = 0, \end{matrix}

(A1)

where t is time and (w_x,w_z) are the components of the acoustic particle velocity. The source term is given by

Q (x, z, t) = \frac{q_{0}}{π B^{2}} exp [- \frac{x^{2} + {(z - z_{s})}^{2}}{B^{2}}] sin (ω t),

(A2)

where

q_{0} = 1

m²s⁻¹ is the source flow rate and B = 0.24 m is its radius.

The FDTD solution uses the optimized fourth-order finite-difference technique (see Ref. 29 for details). The mesh grid is uniform with a grid step of 0.1 m, which corresponds to 17 points per wavelength for f = 200 Hz. The Courant-Friedrichs-Lewy number is set to 0.5, yielding a time step $1.47 \times 10^{- 4}$ s. The FDTD simulation runs for 1.54 s.

The root mean square (rms) sound pressure is determined from the FDTD solution using the last three periods of the sound signal. For comparison with the parabolic equations, the rms pressure is normalized using the factor $q_{0} ω ϱ exp {- ω^{2} B^{2} / [4 c^{2} (z_{s})]}$ ⁠.

REFERENCES

1.

E.

Salomons

,

Computational Atmospheric Acoustics

(

Kluwer Academic

,

Dordrecht

,

2001

), p.

335

.

Google Scholar

Crossref

2.

V. E.

Ostashev

and

D. K.

Wilson

,

Acoustics in Moving Inhomogeneous Media

, 2nd ed. (

CRC Press

,

Boca Raton, FL

,

2015

),

521

pp.

Google Scholar

Crossref

3.

J. F.

Lingevitch

,

M. D.

Collins

,

D. K.

Dacol

,

D. P.

Drob

,

J. C. W.

Rogers

, and

W. L.

Siegmann

, “

A wide angle and high Mach number parabolic equation

,”

J. Acoust. Soc. Am.

111

,

729

–

734

(

2002

).

https://doi.org/10.1121/1.1430683

Google Scholar

Crossref

PubMed

4.

C.

Khodr

,

M.

Azarpeyvand

, and

D. N.

Green

, “

An iterative three-dimensional parabolic equation solver for propagation above irregular boundaries

,”

J. Acoust. Soc. Am.

148

(

2

),

1089

–

1100

(

2020

).

https://doi.org/10.1121/10.0001766

Google Scholar

Crossref

PubMed

5.

C. M.

Nyborg

,

K.

Bolin

,

I.

Karasalo

, and

A.

Fischer

, “

An inter-model comparison of parabolic equation methods for sound propagation from wind turbines

,”

J. Acoust. Soc. Am.

154

,

1299

–

1314

(

2023

).

https://doi.org/10.1121/10.0020562

Google Scholar

Crossref

PubMed

6.

F. D.

Tappert

, “

The parabolic approximation method

,” in

Wave Propagation in Underwater Acoustics.

Lecture Notes in Physics No. 70, edited by

J. B.

Keller

and

J. S.

Papadakis

(

Springer

,

New York

,

1977

).

Google Scholar

Crossref

7.

F. B.

Jensen

,

W. A.

Kuperman

,

M. B.

Porter

, and

H.

Schmidt

,

Computational Ocean Acoustics

, 2nd ed. (

Springer

,

New York

,

2011

), pp.

457

–

529

.

Google Scholar

Crossref

8.

M. D.

Collins

, “

A split-step Padé solution for the parabolic equation method

,”

J. Acoust. Soc. Am.

93

,

1736

–

1742

(

1993

).

https://doi.org/10.1121/1.406739

Google Scholar

Crossref

9.

F.

Dagrau

,

M.

Rénier

,

R.

Marchiano

, and

F.

Coulouvrat

, “

Acoustic shock wave propagation in a heterogeneous medium: A numerical simulation beyond the parabolic approximation

,”

J. Acoust. Soc. Am.

130

,

20

–

32

(

2011

).

https://doi.org/10.1121/1.3583549

Google Scholar

Crossref

PubMed

10.

L.-J.

Gallin

,

M.

Rénier

,

É.

Gaudard

,

T.

Farges

,

R.

Marchiano

, and

F.

Coulouvrat

, “

One-way approximation for the simulation of weak shock wave propagation in atmospheric flows

,”

J. Acoust. Soc. Am.

135

,

2559

–

2570

(

2014

).

https://doi.org/10.1121/1.4869685

Google Scholar

Crossref

PubMed

11.

P. V.

Yuldashev

,

I. S.

Mezdrokhin

, and

V. A.

Khokhlova

, “

Wide-angle parabolic approximation for modeling high-intensity fields from strongly focused ultrasound transducers

,”

Acoust. Phys.

64

(

3

),

309

–

319

(

2018

).

https://doi.org/10.1134/S1063771018030168

Google Scholar

Crossref

12.

R.

Janaswamy

, “

A curvilinear coordinate-based split-step parabolic equation method for propagation predictions over terrain

,”

IEEE Trans. Antennas Propagat.

46

(

7

),

1089

–

1097

(

1998

).

https://doi.org/10.1109/8.704813

Google Scholar

Crossref

13.

V. E.

Ostashev

,

D.

Juvé

, and

P.

Blanc-Benon

, “

Derivation of a wide-angle parabolic equation for sound waves in inhomogeneous moving media

,”

Acta Acust. Acust.

83

(

3

),

455

–

460

(

1997

).

Google Scholar

14.

P.

Blanc-Benon

,

L.

Dallois

, and

D.

Juvé

, “

Long range sound propagation in a turbulent atmosphere within the parabolic approximation

,”

Acta Acust. Acust.

87

(

6

),

659

–

669

(

2001

).

Google Scholar

15.

V. E.

Ostashev

,

D. K.

Wilson

, and

M. B.

Muhlestein

, “

Wave and extra-wide-angle parabolic equations for sound propagation in a moving atmosphere

,”

J. Acoust. Soc. Am.

147

(

6

),

3969

–

3984

(

2020

).

https://doi.org/10.1121/10.0001397

Google Scholar

Crossref

PubMed

16.

D. K.

Wilson

,

M. J.

Shaw

,

V. E.

Ostashev

,

M. B.

Muhlestein

,

M. E.

Swearingen

, and

S. L.

McComas

, “

Numerical modeling of mesoscale infrasound propagation in the Arctic

,”

J. Acoust. Soc. Am.

151

(

1

),

138

–

157

(

2022

).

https://doi.org/10.1121/10.0009163

Google Scholar

Crossref

PubMed

17.

B.

Kayser

,

D.

Mascarenhas

,

B.

Cotté

,

D.

Ecotiere

, and

B.

Gauvreau

, “

Validity of the effective sound speed approximation in parabolic equation models for wind turbine noise propagation

,”

J. Acoust. Soc. Am.

153

(

3

),

1846

–

1854

(

2023

).

https://doi.org/10.1121/10.0017653

Google Scholar

Crossref

PubMed

18.

D.

Mascarenhas

,

B.

Cotté

, and

O.

Doaré

, “

Propagation effects in the synthesis of wind turbine aerodynamic noise

,”

Acta Acust.

7

(23),

1

–

16

(

2023

).

https://doi.org/10.1051/aacus/2023018

Google Scholar

Crossref

19.

J.

Colas

,

A.

Emmanuelli

,

D.

Dragna

,

P.

Blanc-Benon

,

B.

Cotte

, and

R. J. A. M.

Stevens

, “

Wind turbine sound propagation: Comparison of a linearized Euler equations model with parabolic equation methods

,”

J. Acoust. Soc. Am.

154

(

3

),

1413

–

1426

(

2023

).

https://doi.org/10.1121/10.0020834

Google Scholar

Crossref

PubMed

20.

A. D.

Pierce

, “

Wave equation for sound in fluids with unsteady inhomogeneous flow

,”

J. Acoust. Soc. Am.

87

(

6

),

2292

–

2299

(

1990

).

https://doi.org/10.1121/1.399073

Google Scholar

Crossref

21.

V. E.

Ostashev

,

M. B.

Muhlestein

, and

D. K.

Wilson

, “

Extra-wide-angle parabolic equations in motionless and moving media

,”

J. Acoust. Soc. Am.

145

(

2

),

1031

–

1047

(

2019

).

https://doi.org/10.1121/1.5091011

Google Scholar

Crossref

PubMed

22.

H. S.

Wall

,

Analytic Theory of Continued Fractions

(

Chelsea Publishing

,

New York

,

1973

).

Google Scholar

23.

M. D.

Collins

, “

Applications and time-domain solution of higher-order parabolic equations in underwater acoustics

,”

J. Acoust. Soc. Am.

86

(

3

),

1097

–

1102

(

1989

).

https://doi.org/10.1121/1.398101

Google Scholar

Crossref

24.

L.

Fishman

and

J. J.

McCoy

, “

Derivation and application of extended parabolic wave theories. I. The factorized Helmholtz equation

,”

J. Math. Phys.

25

,

285

–

296

(

1984

).

https://doi.org/10.1063/1.526149

Google Scholar

Crossref

25.

A.

Bamberger

,

B.

Engquist

,

L.

Halpern

, and

P.

Joly

, “

Higher order paraxial wave equation approximations in heterogeneous media

,”

SIAM J. Appl. Math.

48

,

129

–

154

(

1988

).

https://doi.org/10.1137/0148006

Google Scholar

Crossref

26.

F.

Collino

, “

Perfectly matched absorbing layers for the paraxial equations

,”

J. Comput. Phys.

131

(

1

),

164

–

180

(

1997

).

https://doi.org/10.1006/jcph.1996.5594

Google Scholar

Crossref

27.

V. E.

Ostashev

,

D. K.

Wilson

,

L.

Liu

,

D. F.

Aldridge

,

N. P.

Symons

, and

D. H.

Marlin

, “

Equations for finite-difference, time-domain simulation of sound propagation in moving inhomogeneous media and numerical implementation

,”

J. Acoust. Soc. Am.

117

(

2

),

503

–

517

(

2005

).

https://doi.org/10.1121/1.1841531

Google Scholar

Crossref

PubMed

28.

J.

Assink

,

R.

Waxler

, and

D.

Velea

, “

A wide-angle high Mach number modal expansion for infrasound propagation

,”

J. Acoust. Soc. Am.

141

(

3

),

1781

–

1792

(

2017

).

https://doi.org/10.1121/1.4977578

Google Scholar

Crossref

PubMed

29.

D.

Dragna

and

P.

Blanc-Benon

, “

Towards realistic simulations of sound radiation by moving sources in outdoor environments

,”

Int. J. Aeroacoust.

13

(

5-6

),

405

–

426

(

2014

).

https://doi.org/10.1260/1475-472X.13.5-6.405

Google Scholar

Crossref

2024

U.S. Government

Phase-preserving narrow- and wide-angle parabolic equations for sound propagation in moving media^a)

I. INTRODUCTION

II. PARABOLIC EQUATIONS

A. EWAPE

B. NAPE

C. WAPE

III. GEOMETRICAL ACOUSTICS

A. Dispersion relation and eikonal equation

B. EWAPE

C. NAPE

D. WAPE

IV. STRATIFIED MOVING MEDIUM

A. Exact equation

B. EWAPE

C. NAPE

D. WAPE

V. NUMERICAL IMPLEMENTATION OF THE PARABOLIC EQUATIONS

A. NAPE

B. WAPE

C. Intercomparison of the parabolic equations

VI. NUMERICAL RESULTS

A. Uniformly moving medium

B. Linear sound speed profile

C. Stratified moving medium

VII. CONCLUSIONS

ACKNOWLEDGMENTS

AUTHOR DECLARATION

Conflict of Interest

DATA AVAILABILITY

Appendix

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

Phase-preserving narrow- and wide-angle parabolic equations for sound propagation in moving mediaa) Free

I. INTRODUCTION

II. PARABOLIC EQUATIONS

A. EWAPE

B. NAPE

C. WAPE

III. GEOMETRICAL ACOUSTICS

A. Dispersion relation and eikonal equation

B. EWAPE

C. NAPE

D. WAPE

IV. STRATIFIED MOVING MEDIUM

A. Exact equation

B. EWAPE

C. NAPE

D. WAPE

V. NUMERICAL IMPLEMENTATION OF THE PARABOLIC EQUATIONS

A. NAPE

B. WAPE

C. Intercomparison of the parabolic equations

VI. NUMERICAL RESULTS

A. Uniformly moving medium

B. Linear sound speed profile

C. Stratified moving medium

VII. CONCLUSIONS

ACKNOWLEDGMENTS

AUTHOR DECLARATION

Conflict of Interest

DATA AVAILABILITY

Appendix

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

Phase-preserving narrow- and wide-angle parabolic equations for sound propagation in moving media^a)