Acoustic Wave Model

Ewert and Schröder presented a family of Acoustic Perturbation Equations (APE) for the simulation of flow-induced acoustic fields. These equations describe noise source terms that are determined by unsteady simulations of compressible or incompressible problems. As Lighthill derived an acoustic analogy from the basic equations of Navier-Stokes, an acoustic wave equation can be derived starting from the acoustic perturbation equations as well.

This model calculates the acoustic pressure history on the whole domain, including the effects of convection and refraction.

See [887] and [901].

The complete system of the acoustic perturbation equations in terms of the perturbation variables $(p^{'}, u^{a})$ is given as

Figure 1. EQUATION_DISPLAY

\frac{\partial p^{'}}{\partial t} + c^{2} \nabla\cdot (\bar{ρ} u^{a} + \bar{v} \frac{p^{'}}{c^{2}}) \approx 0

(4711)

Figure 2. EQUATION_DISPLAY

\frac{\partial u^{a}}{\partial t} + \nabla (\bar{v} \cdot u^{a}) + \nabla (\frac{p^{'}}{\bar{ρ}}) \approx \nabla Φ_{p}

(4712)

where:

$p^{'}$ is the perturbation pressure
$u^{a}$ is the irrotational perturbation velocity
$\bar{ρ}$ is the time-averaged density
$\bar{v}$ is the time-averaged (mean) velocity
$c$ is the speed of sound
$Φ_{p}$ is the noise source function

Using the relationship $p^{'} = \bar{ρ} Φ_{p} + p^{a}$ (where $p^{a}$ is the acoustic pressure) from Ewert and Schroder and taking the substantial derivative $(\frac{\partial}{\partial t} + \bar{v} \cdot \nabla)$ of Eqn. (4711) and the divergence of Eqn. (4712), you can obtain two equations in terms of $p^{a}$ and $u^{a}$ as dependent variables. For incompressible flows, (that is, $\nabla\cdot \bar{v} = 0$ ) and combining the two equations, a single equation is obtained in terms of acoustic pressure and noise sources from incompressible flow:

Figure 3. EQUATION_DISPLAY

\frac{1}{c^{2}} \frac{\partial^{2} p^{a}}{\partial t^{2}} + \frac{2 \bar{v}}{c^{2}} \cdot \nabla \frac{\partial p^{a}}{\partial t} + \frac{(\bar{v} \cdot \nabla)}{c^{2}} (\nabla \cdot \bar{v} p^{a}) - \nabla^{2} p^{a} = - (\frac{\bar{ρ}}{c^{2}} \frac{\partial^{2} Φ_{p}}{\partial t^{2}} + \frac{2 (\bar{v} \cdot \nabla)}{c^{2}} \bar{ρ} \frac{\partial Φ_{p}}{\partial t} + \frac{(\bar{v} \cdot \nabla)}{c^{2}} (\nabla \cdot \bar{v} Φ_{p}))

(4713)

An equation for a wave caused mainly by convection effects is derived from Eqn. (4713) by assuming that, for noise sources from incompressible flows [887], $\bar{ρ} \frac{\partial Φ_{p}}{\partial t} \approx \frac{\partial P^{'}}{\partial t}$ and $\nabla Φ_{p} \approx \frac{1}{\bar{ρ}} \nabla P^{'}$ , where $P^{'}$ is the hydrodynamic pressure fluctuation. In this equation, a physical damping mechanism is necessary to eliminate spurious waves that originate from mesh-coarsening transitions. From the concept of wave propagation in a viscous medium, the final equation with the physical damping term is:

Figure 4. EQUATION_DISPLAY

\frac{1}{c^{2}} \frac{\partial^{2} p^{a}}{\partial t^{2}} + \frac{2 \bar{v}}{c^{2}} \cdot \nabla \frac{\partial p^{a}}{\partial t} + \frac{(\bar{v} \cdot \nabla)}{c^{2}} (\nabla \cdot \bar{v} p^{a}) - \nabla^{2} (p^{a} + τ \frac{\partial p^{a}}{\partial t}) = - (\frac{1}{c^{2}} \frac{\partial^{2} P^{'}}{\partial t^{2}} + \frac{2 (\bar{v} \cdot \nabla)}{c^{2}} \frac{\partial P^{'}}{\partial t} + \frac{(\bar{v} \cdot \nabla)}{c^{2}} (\nabla \cdot \bar{v} P^{'}))

(4714)

where, in the physical damping term,

τ = χ \frac{Δ t}{π λ}

. In this term,

$χ$ is the damping coefficient (= 0 for no damping and = 1 for maximum damping)
$Δ t$ is the time-step
$λ = \frac{c Δ t}{Δ x}$ is the local CFL number (where $Δ x$ is the grid spacing)

Acoustic waves are damped in poorly resolved zones (coarse mesh) using this damping coefficient to suppress undesired oscillations that pollute predictions inside the region of interest.

In the absence of mean flow, Eqn. (4714) reduces to

Figure 5. EQUATION_DISPLAY

\frac{1}{c^{2}} \frac{\partial^{2} p^{a}}{\partial t^{2}} - \nabla^{2} (p^{a} + τ \frac{\partial p^{a}}{\partial t}) = - \frac{1}{c^{2}} \frac{\partial^{2} P^{'}}{\partial t^{2}}

(4715)

This equation resembles Lighthill's equation except that this one is based on noise sources coming from acoustic perturbation equations for incompressible flow.

Low-Pass Time Filters

The sources of noise in the acoustic wave model come mainly from the second derivative of the flow pressure fluctuations with respect to time. Due to mesh size transitions, these sources include high-frequency spurious numerical signals which, in turn, generate high-frequency noise. A suitable low-pass filter suppresses this undesired noise.

For the following filters, the higher the order of the filter, the more signal damping it provides, while also introducing more phase error (that is, the shift in the angle of the wave from that which is expected). You can consider phase error as a measure of the distortion in the signal that is introduced by phase differences for different frequencies.

The Acoustic Wave model applies the low-pass time filter directly to the APE sources. The filtered APE source is:

- \tilde{\frac{1}{c^{2}} \frac{\partial^{2}}{\partial t^{2}} P^{'}} = \frac{1}{c^{2}} \frac{\partial^{2} p^{a}}{\partial t^{2}} - \nabla \cdot \nabla (p^{a} + τ \frac{\partial p^{a}}{\partial t})

First-Order IIR Low-Pass Time Filter

The filtered form of

\frac{d^{2} P^{'}}{d t^{2}}

at time-step

n

is:

Figure 6. EQUATION_DISPLAY

{\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n} = α {(\frac{d^{2} P^{'}}{d t^{2}})}^{n} + (1 - α) {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 1}

(4716)

where:

α = \frac{2 π Δ t f_{c u t}}{2 π Δ t f_{c u t} + 1}

so that the filtered value $\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}$ at a given time-step depends on the filtered value at the previous time-step.

Second-Order Butterworth Low-Pass Time Filter

The filtered form of

\frac{d^{2} P^{'}}{d t^{2}}

at time-step

n

is:

Figure 7. EQUATION_DISPLAY

{\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n} = [\frac{n_{0}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n} + \frac{n_{1}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n - 1} + \frac{n_{2}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n - 2}] + [- \frac{d_{1}}{d_{0}} {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 1} - \frac{d_{2}}{d_{0}} {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 2}]

(4717)

where:

$n_{0} = 1$ , $n_{1} = 2$ , and $n_{2} = 1$
$d_{0} = e^{2} + \sqrt{2} e + 1$
$d_{1} = - 2 (e^{2} - 1)$
$d_{2} = e^{2} - \sqrt{2} e + 1$
$e = \frac{1}{\tan (\frac{2 π Δ t f_{c u t}}{2})}$

Third-Order Butterworth Low-Pass Time Filter

The filtered form of

\frac{d^{2} P^{'}}{d t^{2}}

at time-step

n

is:

Figure 8. EQUATION_DISPLAY

{\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n} = [\frac{n_{0}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n} + \frac{n_{1}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n - 1} + \frac{n_{2}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n - 2} + \frac{n_{3}}{d_{0}} {(\frac{d^{2} P^{'}}{d t^{2}})}^{n - 3}] + [- \frac{d_{1}}{d_{0}} {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 1} - \frac{d_{2}}{d_{0}} {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 2} - \frac{d_{3}}{d_{0}} {\tilde{(\frac{d^{2} P^{'}}{d t^{2}})}}^{n - 3}]

(4718)

where:

$n_{0} = 1$ , $n_{1} = 3$ , $n_{2} = 3$ , and $n_{3} = 1$
$d_{0} = e^{3} + 2 e^{2} + 2 e + 1$
$d_{1} = - 3 e^{3} - 2 e^{2} + 2 e + 3$
$d_{2} = 3 e^{3} - 2 e^{2} - 2 e + 3$
$d_{3} = - e^{3} + 2 e^{2} - 2 e + 1$
$e = \frac{1}{\tan (\frac{2 π Δ t f_{c u t}}{2})}$

Non-Reflective Boundary Condition

The non-reflective boundary condition in the Acoustic Wave model allows acoustic waves to leave the domain through boundary without any spurious reflections. Assuming no source at the boundary and applying the operator " $\nabla\cdot$ " on the momentum equation Eqn. (4712):

Figure 9. EQUATION_DISPLAY

\begin{array}{l} \nabla\cdot \frac{\partial u^{a}}{\partial t} + (\bar{v} \cdot \nabla) \nabla\cdot u^{a} + \frac{1}{\bar{ρ}} \nabla\cdot \nabla p^{a} = 0 \\ \nabla\cdot u^{a} = - \frac{1}{\bar{ρ} c^{2}} (\frac{\partial p^{a}}{\partial t} + \bar{v} \cdot \nabla p^{a}) \end{array}

(4719)

The acoustic pressure is related to the normal component of the acoustic velocity:

Figure 10. EQUATION_DISPLAY

p^{a} = \bar{ρ} c u^{a} \cdot n

(4720)

Taking the time derivative of Eqn. (4720) and using Eqn. (4719), the gradient of the acoustic pressure then becomes:

Figure 11. EQUATION_DISPLAY

\nabla p^{a} \cdot A_{f} = - \frac{(1 - \bar{v} \cdot n / c) \partial p^{a}}{(1 - {\bar{v}}^{2} / c^{2}) c \partial t} ‖ A_{f} ‖

(4721)

where $A_{f}$ is the face normal area and $‖ A_{f} ‖$ is its magnitude.

Newmark Alpha Method

The Newmark Alpha method starts with the following expression for time derivatives of acoustic pressure:

M {(\frac{\partial^{2} p^{a}}{\partial t^{2}})}^{n + 1} + C {(\frac{\partial^{} p^{a}}{\partial t^{}})}^{n + 1} + K {(p^{a})}^{n + 1} = F^{n + 1}

where M, C, and K are operators:

$M = \frac{1}{c^{2}} \int_{V} d V$ is the operator for mass.
$C = \frac{1}{c^{2}} \int_{s} \bar{v} d s$ is the operator for damping.
$K = \int_{s} \nabla\cdot d s + \frac{1}{c^{2}} (\bar{v} \cdot \nabla) \bar{v} \cdot d s$ is the operator for stiffness.
$F$ is the noise source term.
$n$ is the previous time-step.
$n + 1$ is the time-step for which the solution is being computed.

The method then discretizes the above equation and introduces the parameter $α$ to obtain the following modified governing equation:

Figure 12. EQUATION_DISPLAY

M {(\frac{\partial^{2} p^{a}}{\partial t^{2}})}^{n + 1} + (1 + α) C {(\frac{\partial^{} p^{a}}{\partial t^{}})}^{n + 1} + (1 + α) K {(p^{a})}^{n + 1} = {(1 + α) F}^{n + 1} - α F^{n} + α C {(\frac{\partial^{} p^{a}}{\partial t^{}})}^{n} + α K {(p^{a})}^{n}

(4722)

The operators and the external force take the form of matrices, where:

$M$ is the mass matrix
$C$ is the damping matrix
$K$ is the stiffness matrix
$F$ is the external force matrix
$α$ is a weighting parameter, $- 0.33 \leq α \leq 0$ .

APE Fields for the Post FW-H Model

The Acoustic Wave model can be combined with the Ffowcs Williams-Hawkings model to calculate the following fields for far-field cases:

Acoustic Density $ρ^{a}$ :
Figure 13. EQUATION_DISPLAY

$ρ^{a} = \frac{p^{a}}{c^{2}} + ρ_{f}$

(4723)
Acoustic Velocity $u^{a}$ :
$u^{a}$ at the nth time-step, $u_{n}^{a}$ , is computed from the acoustic velocities of the previous time-step and the acoustic pressure gradient. This is obtained from second-order backward differencing of the momentum equation of the Acoustic Perturbation Equations.
Figure 14. EQUATION_DISPLAY

$u_{n}^{a} = \frac{1}{3} [4 u_{n - 1}^{a} - u_{n - 2}^{a} - \frac{2 Δ t \nabla p^{a}}{ρ_{0}}]$

(4724)
APE Total Density $Σ ρ$ :
Figure 15. EQUATION_DISPLAY

$Σ ρ = \frac{P^{a}}{c^{2}} + ρ_{f}$

(4725)
APE Total Velocity $Σ v$ :
Figure 16. EQUATION_DISPLAY

$Σ v = v^{a} + v_{f}$

(4726)
APE Total Pressure $Σ P$ :
Figure 17. EQUATION_DISPLAY

$Σ P = P^{a} + P_{f}$

(4727)

where:

$c$ is the speed of sound.
$ρ_{f}$ is the density of the incompressible flow.
$P^{a}$ is the acoustic pressure.
$P_{f}$ is the pressure of the incompressible flow.
$v_{f}$ is the velocity of the incompressible flow.

Smoothing Algorithm for the Acoustic Pressure Reconstruction Gradient

The reconstructed gradient of acoustic pressure is used to calculate the acoustic pressure at the face centroids. This reconstructed gradient can be smoothed out, eliminating local spikes, to assist convergence. The smoothing algorithm uses the following steps:

Compute the locally averaged normalized and smoothed gradient $\bar{\nabla (p_{c}^{a})}$ :
Figure 18. EQUATION_DISPLAY

$\bar{\nabla (p_{c}^{a})} = \frac{\frac{\nabla (p_{c}^{a})}{{‖ \nabla (p_{c}^{a}) ‖}^{2}} + \sum_{k \in {N B}_{c}}^{} \frac{\nabla (p_{k}^{a})}{{‖ \nabla (p_{k}^{a}) ‖}^{2}}}{\frac{1}{{‖ \nabla (p_{c}^{a}) ‖}^{2}} + \sum_{k \in {N B}_{c}}^{} \frac{1}{{‖ \nabla (p_{k}^{a}) ‖}^{2}}}$

(4728)
where ${N B}_{c}$ is the set of those cells that are direct neighbors, through faces, of cell $c$ .
Copy the smoothed gradient $\bar{\nabla (p_{c}^{a})}$ into the working gradient $\nabla (p_{c}^{a})$ .
Repeat steps 1 and 2 as many times as needed to eliminate local spikes.