**3. Numerical procedure**

#### **3.1 Transient CFD simulation**

Transient flow fields around rod-airfoil and airfoil-airfoil models are simulated at Reynolds number *Red* = 28,800 and *Rec* = 288,000; those are based on the rod diameter and the airfoil chord, respectively. The inflow velocity in *x* direction, temperature, and density are *U*<sup>∞</sup> *=* 72 m/s,*T*<sup>∞</sup> *=* 300 K, and *ρ*<sup>∞</sup> *=* 1.177 kg/m3 , respectively and Mach number is *M* = 0.207. **Figure 3** shows a computational domain for the rod-airfoil and the airfoil-airfoil models. For these simulations, the three-dimensional computational domain has been applied, as shown in **Figure 3**. Unsteady flow fields are calculated using the commercial CFD code ANSYS Fluent version 2019R1 and its compressible LES (Dynamic Smagorinsky model) calculation features. Steady velocities are imposed on the inflow boundary. Pressure boundary conditions are applied on the top, bottom, and outflow boundaries. Nonreflecting boundary conditions are applied on top, bottom, inflow, and outflow boundaries. No-slip conditions are applied on the walls. At the boundaries at *z =* 0.009 m and 0.009 m, the periodic conditions are applied.

The domains of the rod-airfoil simulation model contain 1,169,322 and 1,186,372 hex cells for *L/d =* 0.2 and 1.0, respectively. The domains of the airfoil-airfoil simulation model contain 1,423,114; 1,435,343; and 1,434,356 hex cells for *L/c =* 0.2,

**Figure 3.** *Computational domain.*

*Wake-Body Interaction Noise Simulated by the Coupling Method Using CFD and BEM DOI: http://dx.doi.org/10.5772/intechopen.92783*

**Figure 4.** *Computational mesh. (a) Rod-airfoil model; (b) airfoil-airfoil model.*

0.6, and 1.0, respectively. **Figure 4(a)** and **(b)** shows the computational meshes near the rod-airfoil and the airfoil-airfoil models, respectively. The cell spacing adjacent to the wall is 0.2 mm (0.033*d*). Steady-state simulations were performed using Spalart-Allmaras (S-A) turbulence model and then used as initial conditions of transient LESs. The transient simulations were performed for 50,000 time steps with a time step size *Δt =* 2e-6 s.

#### **3.2 Lighthill equation**

prescribed to be 0.2, 0.6, and 1.0, and those are equal to *L/d =* 2.0, 6.0, and 10. In this chapter, the Reynolds number is fixed to be *Rec* = 288,000, which is based on the airfoil chord *c*. The spanwise length of the airfoil and the airfoil is 0.3*c,* which

*Schematic diagram of arifoil-airfoil model. (a) Airfoil-airfoil model; (b) parameters.*

Transient flow fields around rod-airfoil and airfoil-airfoil models are simulated at Reynolds number *Red* = 28,800 and *Rec* = 288,000; those are based on the rod diameter and the airfoil chord, respectively. The inflow velocity in *x* direction, temperature, and density are *U*<sup>∞</sup> *=* 72 m/s,*T*<sup>∞</sup> *=* 300 K, and *ρ*<sup>∞</sup> *=* 1.177 kg/m3

The domains of the rod-airfoil simulation model contain 1,169,322 and 1,186,372

hex cells for *L/d =* 0.2 and 1.0, respectively. The domains of the airfoil-airfoil simulation model contain 1,423,114; 1,435,343; and 1,434,356 hex cells for *L/c =* 0.2,

respectively and Mach number is *M* = 0.207. **Figure 3** shows a computational domain for the rod-airfoil and the airfoil-airfoil models. For these simulations, the three-dimensional computational domain has been applied, as shown in **Figure 3**. Unsteady flow fields are calculated using the commercial CFD code ANSYS Fluent version 2019R1 and its compressible LES (Dynamic Smagorinsky model) calculation features. Steady velocities are imposed on the inflow boundary. Pressure boundary conditions are applied on the top, bottom, and outflow boundaries. Nonreflecting boundary conditions are applied on top, bottom, inflow, and outflow boundaries. No-slip conditions are applied on the walls. At the boundaries at

*z =* 0.009 m and 0.009 m, the periodic conditions are applied.

,

corresponds to 3*d*.

**Figure 2.**

**Figure 3.**

**124**

*Computational domain.*

**3. Numerical procedure**

**3.1 Transient CFD simulation**

*Vortex Dynamics Theories and Applications*

Lighthill equation [24, 25] in the frequency domain is derived from the equation of continuity and compressible Navier-Stokes equation and as follows:

$$(\nabla^2 + k^2)p = -\frac{\partial^2 T\_{lm}}{\partial \mathbf{x}\_l \partial \mathbf{x}\_m} \tag{1}$$

where *p* is the acoustic pressure, *k* is the wave number, *c*<sup>∞</sup> is the speed of sound, *l* and *m* indicate each direction in the Cartesian coordinates, and *v* is the flow velocity. *Tlm* is the Lighthill stress tensor and as follows:

$$T\_{lm} = \rho v\_l \upsilon\_m + (p - c\_\infty \,^2 \rho) \delta\_{lm} - \tau\_{lm},\tag{2}$$

where ρ is the density and is 1.225 kg/m3, *δij* is the Kronecker delta, and *τlm* is the viscous stress tensor. For a low-Mach number and high-Reynolds number flow regime, the second and third terms of Eq. (2) are negligible [26–29]. Therefore, the first term is used for the present work.

#### **3.3 Extraction of acoustic source**

To convert the acoustic source time histories into the frequency spectra, the discrete Fourier transform (DFT) has been applied. The acoustics sources are extracted from 1250 steps (from *t =* 0.05 s to 0.1 s) flow field data, the sampling period is 4e-5 s.

#### **3.4 Acoustic simulation**

The BEM solver in commercial acoustic simulation package, WAON, is used to solve the acoustic characteristics [30]. In a sound field that satisfies the threedimensional Helmholtz equation, the Kirchhoff-Helmholtz integral equation [31]

for sound pressure is described as follows with respect to a point i and an area S of a surface on a boundary.

$$\frac{1}{2}P(r\_i) = \int\_{\varGamma} \left( P(r\_q) \frac{\partial G(r\_i, r\_q)}{\partial n\_q} - \frac{\partial P(r\_q)}{\partial n\_q} G(r\_i, r\_q) \right) d\mathbf{S} + p\_d(r\_q) \tag{3}$$

In this solver, the following simultaneous linear equation is solved:

$$(\mathbf{E} + \mathbf{B} + \mathbf{C})\mathbf{p} = j\alpha\rho\mathbf{A} + \mathbf{p}\_d \tag{4}$$

Here, p is the acoustic pressure vector, *v* is the particle velocity vector, and the entries of the influence coefficient matrices are represented as follows:

$$E\_{\vec{\eta}} = \frac{1}{2} \delta\_{\vec{\eta}},\tag{5}$$

**Figure 5.**

**Figure 6.**

**127**

*present (*Red = *28,800).*

*Boundary elements. (a) Rod-airfoil model; (b) airfoil-airfoil model.*

*DOI: http://dx.doi.org/10.5772/intechopen.92783*

*Wake-Body Interaction Noise Simulated by the Coupling Method Using CFD and BEM*

*Mean velocity and RMS of fluctuation velocity in* x *direction distributions at* x/c = *0.87 (a), (b) and 0.25 (c), (d). (a) Velocity; (b) RMS value of fluctuation velocity. Measured by Jacob et al. [9]*

*(*Red = *48,000); Agawal and Sharma [18] (*Red = *48,000); Jiang et al. [21] (*Red = *48,000); present (*Red = *28,800). (c) Velocity; (d) RMS value of fluctuation velocity. Measured by Jacob et al. [9] (*Red = *48,000); Boudet et al. [16] (*Red = *48,000); Jiang et al. [21] (*Red = *48,000);*

$$A\_{\vec{\eta}} = \int\_{\Gamma\_1} N\_j(\mathbf{r}\_q) G(\mathbf{r}\_i, \mathbf{r}\_q) d\mathbf{S}\_q,\tag{6}$$

$$B\_{\vec{\eta}} = \int\_{\varGamma} N\_{\vec{\jmath}}(\mathbf{r}\_q) \frac{\partial G(\mathbf{r}\_i, \mathbf{r}\_q)}{\partial \mathbf{n}\_q} dS\_q,\tag{7}$$

$$\mathbf{C}\_{ij} = \frac{jk}{Z\_j} \int\_{\Gamma\_2} N\_j(\mathbf{r}\_q) G(\mathbf{r}\_i, \mathbf{r}\_q) d\mathbf{S}\_q,\tag{8}$$

$$G(r\_i, r\_q) = \frac{\mathbf{e}^{jk \, |r\_p - r\_i|}}{4\pi |r\_p - r\_s|} \tag{9}$$

where *δij* is Kronecker delta, and *Γ<sup>1</sup>* is a vibration boundary and a part of *Γ*. *Γ* is the total boundary. *Γ<sup>2</sup>* is an absorption boundary and a part of *Γ*. *ri* is the position vector at the node *i*, **rq** is the position vector of the source point *q,* and *Nj* is the interpolation function of the node *j*. *zj* is the acoustic impedance ratio at the node *j*. *G* is the fundamental solution of a three-dimensional sound field. With the number of nodes *N*, the component *p* of the vector *p* is expressed as follows:

$$p(\mathbf{r}\_q) = \sum\_{j=1}^{N} N\_j(\mathbf{r}\_q) p\_j. \tag{10}$$

The component *pd* of the vector *pd* is the direct pressure contribution from the acoustic source, which is evaluated by the following equation:

$$p\_{\rm d}(r\_{\rm p}) = \frac{1}{4\pi} \frac{\partial^2}{\partial \mathbf{x}\_l \partial \mathbf{x}\_m} \int\_{-\infty}^{\infty} \frac{T\_{lm}(r\_s, \alpha) \mathbf{e}^{jk|r\_p - r\_l|}}{|r\_p - r\_s|} dV \tag{11}$$

where *∂*<sup>2</sup> *=∂xl∂xm* is the directional derivative and *V* is the volume of the flow field (in this case, the region filled by CFD cell). *rp* is a position vector of the monitor point p, and *rs* is a position vector of the source point s. There are 8340 boundary elements. The acoustic sources are extracted from CFD results, whose numbers are equivalent to the number of grids of the CFD model. **Figure 5** shows the boundary elements.

#### **3.5 Validation of numerical results**

The numerical accuracy of the computation has been examined for the case of *L/d =* 10 in the rod-airfoil model. In order to validate the present computation, we *Wake-Body Interaction Noise Simulated by the Coupling Method Using CFD and BEM DOI: http://dx.doi.org/10.5772/intechopen.92783*

**Figure 5.** *Boundary elements. (a) Rod-airfoil model; (b) airfoil-airfoil model.*

#### **Figure 6.**

for sound pressure is described as follows with respect to a point i and an area S of a

� � � �

Here, p is the acoustic pressure vector, *v* is the particle velocity vector, and the

� �*G ri*, *r<sup>q</sup>*

� � *<sup>∂</sup><sup>G</sup> <sup>r</sup>i*, *<sup>r</sup><sup>q</sup>* � � *∂nq*

� � <sup>¼</sup> <sup>e</sup>*jk*j j *<sup>r</sup>p*�*r<sup>s</sup>*

where *δij* is Kronecker delta, and *Γ<sup>1</sup>* is a vibration boundary and a part of *Γ*. *Γ* is the total boundary. *Γ<sup>2</sup>* is an absorption boundary and a part of *Γ*. *ri* is the position vector at the node *i*, **rq** is the position vector of the source point *q,* and *Nj* is the interpolation function of the node *j*. *zj* is the acoustic impedance ratio at the node *j*. *G* is the fundamental solution of a three-dimensional sound field. With the number

*j*¼1

ð<sup>∞</sup> �∞

field (in this case, the region filled by CFD cell). *rp* is a position vector of the monitor point p, and *rs* is a position vector of the source point s. There are 8340 boundary elements. The acoustic sources are extracted from CFD results, whose numbers are equivalent to the number of grids of the CFD model. **Figure 5** shows

The numerical accuracy of the computation has been examined for the case of *L/d =* 10 in the rod-airfoil model. In order to validate the present computation, we

The component *pd* of the vector *pd* is the direct pressure contribution from the

*N <sup>j</sup> r<sup>q</sup>* � �*<sup>p</sup> <sup>j</sup>*

*=∂xl∂xm* is the directional derivative and *V* is the volume of the flow

*Tlm <sup>r</sup><sup>s</sup>* ð Þ ,*<sup>ω</sup>* <sup>e</sup>*jk*j j *<sup>r</sup>p*�*r<sup>s</sup> r<sup>p</sup>* � *r<sup>s</sup>* � � � �

� �*G ri*, *r<sup>q</sup>*

4π *r<sup>p</sup>* � *r<sup>s</sup>* � � � �

*Eij* <sup>¼</sup> <sup>1</sup> 2

� *<sup>∂</sup><sup>P</sup> <sup>r</sup><sup>q</sup>* � � *∂nq*

*G ri*, *r<sup>q</sup>*

ð Þ *E* þ *B* þ *C p* ¼ *jωρA* þ *pd* (4)

*dS* þ *pd rq*

*δij*, (5)

� �*dSq*, (6)

� �*dSq*, (8)

*dSq*, (7)

*:* (10)

*dV* (11)

(9)

� � (3)

surface on a boundary.

1 2

*P*ð Þ¼ *r<sup>i</sup>*

*Vortex Dynamics Theories and Applications*

ð *Γ*

*P r<sup>q</sup>*

� � *<sup>∂</sup><sup>G</sup> <sup>r</sup>i*, *<sup>r</sup><sup>q</sup>* � � *∂nq*

In this solver, the following simultaneous linear equation is solved:

entries of the influence coefficient matrices are represented as follows:

*Aij* ¼ ð Γ1 *N <sup>j</sup> r<sup>q</sup>*

*Bij* ¼ ð *Γ N <sup>j</sup> r<sup>q</sup>*

*Cij* <sup>¼</sup> *jk Z j* ð Γ2 *N <sup>j</sup> r<sup>q</sup>*

*G ri*, *r<sup>q</sup>*

of nodes *N*, the component *p* of the vector *p* is expressed as follows:

*p r<sup>q</sup>* � � <sup>¼</sup> <sup>X</sup> *N*

acoustic source, which is evaluated by the following equation:

*∂*2 *∂xl∂xm*

*p*<sup>d</sup> *r*<sup>p</sup> � � <sup>¼</sup> <sup>1</sup> 4π

where *∂*<sup>2</sup>

**126**

the boundary elements.

**3.5 Validation of numerical results**

*Mean velocity and RMS of fluctuation velocity in* x *direction distributions at* x/c = *0.87 (a), (b) and 0.25 (c), (d). (a) Velocity; (b) RMS value of fluctuation velocity. Measured by Jacob et al. [9] (*Red = *48,000); Agawal and Sharma [18] (*Red = *48,000); Jiang et al. [21] (*Red = *48,000); present (*Red = *28,800). (c) Velocity; (d) RMS value of fluctuation velocity. Measured by Jacob et al. [9] (*Red = *48,000); Boudet et al. [16] (*Red = *48,000); Jiang et al. [21] (*Red = *48,000); present (*Red = *28,800).*

compared the present results with the experimental measurements performed by Jacob et al. [9] and the numerical results performed by Boudet et al. [16], Agawal and Sharma [18], and Jiang et al. [21], respectively. However, it is specified that the present study is performed at *Red =* 28,800, *d =* 0.06 m, and other experiments and the simulations are performed at *Red =* 48,000, *d =* 0.1 m. **Figure 6** shows the mean velocity and RMS value of the fluctuation velocity in *x* (streamwise) direction normalized by the incoming velocity at two locations *x/c =* 0.87 and *x/c =* 0.25. As shown in **Figure 6(a)**, the present calculation predicted similar mean velocity profile in the streamwise direction compared with other numerical results. As shown in **Figure 6(b)**, the RMS value of the fluctuation velocity in the streamwise direction obtained by the present calculation is close to those obtained by the other numerical calculations; however, the difference of the fluctuation velocity near the center-line can be seen. A possible cause for this result is the difference of the methods, meshes and Reynolds number. **Figure 6(c)** and **(d)** shows that the calculated profiles of the mean velocity and the RMS of the fluctuation velocity in the streamwise direction represent a good agreement with those from the experiments and other numerical results.

**4. Results**

**Figure 8.**

**Figure 9.**

**129**

*(d)* L/d = *10.*

*Vorticity in the z direction. (a)* L/d = *2; (b)* L/d = *10.*

*4.1.1 Flow patterns*

**4.1 Rod-airfoil simulation results**

*DOI: http://dx.doi.org/10.5772/intechopen.92783*

Typical examples of instantaneous vorticity fields are presented in **Figure 8**. In the present calculation, for the case of *L/d =* 2, when the spacing between the rod and the airfoil is small, the boundary layers separated from the rod upstream did not roll up and reattaches to the airfoil downstream. The shear layers rolled up and formed vortices on the downstream airfoil surface, and the vortices were shed and convected downstream. The Kármán-vortex street formed from the rod was suppressed for the short spacing, and this mode is called "non-shedding mode" and "the Kármán-vortex street suppressing mode" as indicated by Munetaka et al. [10] and Jiang et al. [21], respectively. On the other hand, for the case of *L/d =* 10, when the spacing between the rod and the airfoil is large, the boundary layers separated from the rod upstream rolled up and formed vortices in the region between the rod and the airfoil, and the formed vortices shed from the rod interacted with the airfoil and impinged on the leading edge of the airfoil. The impinged vortices were distorted and convected downstream. This mode is called "the Kármán-street shedding mode" as indicated by Jiang et al. [21], and also called wake-body interaction

*Wake-Body Interaction Noise Simulated by the Coupling Method Using CFD and BEM*

or body-vortex interaction [6, 16, 21]. Similar phenomena such as the

*Time-averaged velocity in the streamwise (*x*) and the vertical (*y*) directions. (a), (b) velocity in the streamwise (*x*) direction; (c), (d) velocity in the vertical (*y*) direction. (a)* L/d = *2; (b)* L/d = *10; (c)* L/d = *2;*

**Figure 7** shows the spectra of the SPL at the location (*x =* 0.68 m, *y =* 1.74 m) calculated by the acoustic BEM simulation using the acoustic sources extracted from the CFD results. The peak frequency (*St* 0.2) and the spectrum around the peak frequency are well predicted compared with the experimental result obtained by Jacob et al. [9]. The peak value of the SPL obtained by the present study is slightly lower than that from the experimental result. There is a difference in the spectrum at high frequencies, *St >* 0.6. A possible cause for these differences of the peak SPL and the spectrum at high frequencies is that the Reynolds number is different between the present study and the experiments, and the mesh and time resolutions for the calculation might not be enough for the accurate prediction of the spectrum at high frequencies. However, the dominant peak SPL and SPLs around the peak frequency are well predicted in the present calculation. The peak frequency *St* 0.2 as shown in **Figure 7** almost corresponds to the vortex shedding frequency of the cylinder [2], which means that the vortex shedding from the rod plays an important role for the noise generation from the rod-airfoil model.

**Figure 7.** *SPL spectra at the location* x = *0.68 m,* y = *1.74 m. Measured by Jacob et al. [9] (*Red = *48,000); present (*Red = *28,800).*

*Wake-Body Interaction Noise Simulated by the Coupling Method Using CFD and BEM DOI: http://dx.doi.org/10.5772/intechopen.92783*
