Numerical Investigation of Wind Turbine Wake Characteristics under Different Wind Shear Conditions Using the Actuator Line Model

Chengbo Peng, Guoyi He^*, Kaixu Zhang

School of Aeronautics and Astronautics, Nanchang Hangkong University, Nanchang, China.

DOI: 10.4236/jpee.2025.136008                

Abstract

For the ultra-high-power IEA 22 MW wind turbine, an actuator line model was first developed based on actuator line theory. This model was then integrated into a CFD flow field as a momentum source term for numerical simulation. Simulations were performed under different wind shear conditions to examine how the wind shear exponent affects the wake characteristics and power output of the turbine. The findings of this study provide a scientific foundation for wake modeling of ultra-large wind turbines and optimizing wind farm design.

Keywords

Wind Turbine, Large-Eddy Simulation, Actuator Line Model, Wake Analysis

1. Introduction

As a crucial part of renewable energy, wind power is drawing increased attention amid growing concerns over environmental degradation [1]-[3]. Numerical simulation and experimental measurement are the main methods used in wind turbine studies, especially actuator models that can simulate wind farms [4] [5]. Shen et al. [6] used an actuator line model to numerically study the MEXICO experimental wind turbine. Their results showed that the model could precisely predict the wake expansion, vortex radius, circulation, and the axial and tangential wind speed distribution of the MEXICO wind turbine. Zhang et al. [7] used the LES method to simulate a large wind farm with a regular turbine layout that was aligned with the wind direction, and then analyzed the influence of the vertical staggering distance. Onel and Tuncer [8] studied wake recovery and power generation evaluation for turbines in tandem configurations. Under turbulent inflow, velocity deficits in the wake were mitigated, and the power generation of downstream turbines increased. Tu et al. [9] investigated the aerodynamic effects of yaw control on two tandem wind turbines under uniform inflow conditions. They noted that when the upstream rotor yawed, the downstream rotor compensated for the power loss, but this also increased the cyclic loads experienced by the blades. Parada [10] proposed a method to maximize wind farm power and efficiency, proposing that regularly arranged wind farms were not optimally positioned due to wake losses. Furthermore, it was found that Parada’s approach could improve the average efficiencies by 4.09%. Zhu et al. [11] employed dynamic mode decomposition (DMD) to analyze the vortex system downstream of the rotor, indicating that tower shedding vortices accelerate the disintegration of tip vortices and enhance instability.

Research on the wake characteristics of the IEA 22 MW wind turbine remains limited both domestically and internationally. Investigating its wake behavior and power distribution under different wind shear conditions is crucial not only for enhancing the overall efficiency of wind farms but also for optimizing turbine operation and minimizing equipment wear. This study, using the actuator line method, systematically reveals the wake evolution mechanisms and power distribution of the IEA 22 MW ultra-large wind turbine under complex wind shear conditions, providing key technical insights to support the future development of the wind energy industry.

2. Numerical Method and Simulation Setup

2.1. Governing Equation

The continuity and Navier-Stokes equations for the incompressible fluid are shown below:

$\frac{\partial {\overline{u}}_i}{\partial x_i}=0$ (1)

$\frac{\partial {\overline{u}}_i}{\partial t}+\frac{\partial }{\partial x_j}\left({\overline{u}}_i{\overline{u}}_j\right)=-\frac{1}{\rho }\frac{\partial \overline{p}}{\partial x_i}+\frac{\mu }{\rho }\left(\frac{{\partial }^2{\overline{u}}_i}{\partial x_j\partial x_j}-\frac{\partial {\tau }_{ij}}{\partial x_j}\right)+\frac{f_i}{\rho }$ (2)

where the terms from right to left in the Navier–Stokes equation represent the forces exerted by the rotor of turbine, subgrid scale stress, viscosity stress, pressure, convective, and transient, respectively. In this paper, the commonly used Smagorinsky-Lilly sub-grid scale (SGS) model for LES is applied.

2.2. Actuator Line Model

The body force distributed on the blade is

$f_{2D}=\frac{1}{2}\rho U_{rel}^{2}c\left(C_l{e}_l+C_d{e}_d\right)\cdot F$ (3)

where $\rho$ is the air density, $U_{rel}$ is the relative inflow velocity in front of the rotor, $c$ is the chord length, and $C_l$ is the lift coefficient, $C_d$ is the drag coefficient. $e_l$ and $e_d$ are the unit vectors in the directions of lift and drag, respectively. Taking into account the tip loss correction $F_{tip}$ and root loss correction $F_{hub}$ , the overall blade root-tip loss factor $F$ is calculated using the following equations:

$F=F_{tip}\cdot F_{hub}$ (4)

$F_{tip}=\frac{2}{\pi }\arccos \left[\exp \left(-g\frac{B\left(R-r\right)}{2r\sin \phi }\right)\right]$ (5)

$F_{hub}=\frac{2}{\pi }\arccos \left[\exp \left(-\frac{B\left(r-R_{hub}\right)}{2r\sin \phi }\right)\right]$ (6)

here, $g=\exp \left(-0.125\left(B\lambda -21\right)\right)+0.1$ , $B$ is the number of turbine blades, $R$ is the blade tip radius, $R_{hub}$ is the hub radius, and $\lambda$ is the tip speed ratio (TSR).

The Gaussian projection method is applied to smooth the force distribution and convert line forces into a body force field, effectively preventing numerical discontinuity,

$f_{\epsilon }=f_{2D}\otimes {\eta }_{\epsilon }$ (7)

${\eta }_{\epsilon }=\frac{1}{{\epsilon }^3{\pi }^{3/2}}\exp \left(-{\left(\frac{d}{\epsilon }\right)}^2\right)$ (8)

where $d$ is the distance between the grid and the actuator point, and $\epsilon$ is the Gaussian projection width. To avoid numerical oscillations and enhance computational accuracy, the value of $\epsilon$ is set to $2\Delta$ .

Similar to the actuator line model used for the rotor, the nacelle and tower axes are divided into several actuator points. The body forces calculated from the nacelle and tower models are projected using a Gaussian distribution function to simulate the wakes of the nacelle and tower. When the wind turbine adjusts its yaw, the nacelle axis (AL) rotates by the corresponding angle. The body forces of the nacelle and tower are defined as:

$f_{nac}=-\frac{1}{2}\rho U_{nac}^2{C}_{D,nac}{A}_{nac}$ (9)

$f_{tow}=-\frac{1}{2}\rho U_{tow}^2{C}_{D,tow}{d}_{tow}dh$ (10)

where $U_{nac}$ is the axial inflow velocity of the nacelle, $C_{D,nac}$ is the drag coefficient of the nacelle, and $A_{nac}$ is the projected area of the nacelle within its plane of rotation. $U_{tow}$ is the axial inflow velocity of the tower, $C_{D,tow}$ is the drag coefficient of the tower, $d_{tow}$ is the diameter of the circular cross-section at the center of the tower, and $dh$ refers to the height of each element of the tower.

2.3. Simulation Setup

This study utilizes the IEA 22 MW wind turbine, with its detailed key parameters listed in Table 1.

Table 1. Properties of IEA 22 MW wind turbine.

Property	Value
Rated power	22 MW
Hub height	170 m
Rotor diameter	284 m
Hub diameter	8.4 m
Rated wind speed	11.0 m/s
Design tip-speed ratio	9.153
Cut-in wind speed	3 m/s
Cut-out wind speed	25 m/s

Figure 1 illustrates the computational domain. The wind farm scale is based on the diameter ($D$ ) of the IEA 22 MW wind turbine, with the flow direction (x), vertical direction (y), and lateral direction (z) having dimensions of $18D\times 4D\times 5D$ , respectively.

Figure 1. Schematic diagram of computational domain size for IEA 22 MW wind turbines.

All numerical calculations in this study use a Cartesian grid, with grid refinement applied around the two rotor regions. Figure 2 illustrates the wake interaction computational grid. The grid points are evenly distributed in the x, y, and z directions. Outside of the near-field area, the grid spacing increases at a rate no greater than 1.1 in the x, y, and z directions. The time step is set at 0.0236 seconds, which corresponds to 360 time steps for one full rotation of the blades, satisfying the Courant-Friedrichs-Lewy (CFL) condition.

Figure 2. Schematic diagram of computational mesh for IEA 22 MW wind turbine.

This study employs the power-law distribution to model the wind shear effect:

$V\left(y\right)=V_{hub}{\left(\frac{y}{y_{hub}}\right)}^{\alpha }$ (11)

where $V_{hub}$ is the wind speed at hub height, and $y_{hub}$ denotes the hub height. The wind shear exponent $\alpha$ in this study is assigned values of 0.12, 0.15, 0.22, and 0.30. The mean incoming wind speed in the turbulent inflow region is 11.0 m/s, with a turbulence intensity of 6%. In this study, the simulation of the IEA 22 MW wind turbine runs for 764 seconds, during which the turbine completes about 90 rotations.

3. Results

3.1. Main Wake Characteristics

Figure 3 presents instantaneous velocity slices in the vertical plane at the wind turbine location under varying wind shear exponents. It can be seen that the velocity difference across the rotor plane decreases with height, and this effect becomes more pronounced as the wind shear exponent increases. A similar trend is observed in the velocity difference above and below the hub center, suggesting that as the wind shear exponent increases, the wind energy utilization at the upper and lower ends of the rotor becomes more balanced. Additionally, with the influx of external airflow, the wake recovery effect is enhanced, causing the velocity deficit in the wake to gradually disappear around $X=1000\text{m}$ . The larger the wind shear exponent, the earlier this recovery occurs.

Figure 3. Contours of wake velocity.

Figure 4 depicts the vorticity isosurfaces under different wind shear exponents, with wake velocity used for coloring. As the wind shear exponent increases, the spacing between tip vortices also widens. The axial movement of the tip vortices is approximately equal to the local wind speed. Since the top vortices are in a high-speed region, they travel much faster than the bottom vortices, leading to axial stretching of adjacent vortex rings.

Figure 4. Vorticity isosurface colored by wake velocity.

Figure 5 presents vorticity slices in the vertical plane at the wind turbine location under different wind shear exponents. The tip vortex expansion phenomenon can be observed, but it weakens as the wind shear exponent increases. The tip vortices begin to break up around two rotor diameters downstream, with the breakup pattern differing between the top and bottom of the rotor, forming a “longer top, shorter bottom” effect. This occurs because the stronger momentum replenishment in the high-speed flow at the top enhances wake recovery. However, the top vortices are more intense, leading to a deeper velocity deficit and a longer recovery distance. Meanwhile, at the bottom, the momentum difference between the low-speed ambient flow and the wake is smaller, promoting faster turbulence mixing and accelerating wake recovery, causing the vortex structures to break up earlier.

Figure 6 shows the vertical profiles of dimensionless axial velocity at different downstream positions under various wind shear indices. At the near-wake position ($X=1D$ ), the wind turbine velocity profile exhibits asymmetry due to the combined effects of wind shear inflow and tower shadow effects. The velocity profiles under all wind shear indices display a bow-shaped pattern, with smaller velocity deficits observed in the upper rotor height region under lower wind shear indices. At the $X=3D$ position, the wake of the IEA 22 MW wind turbine resides in the wake interaction zone. The velocity deficit within the rotor plane height range shows further intensification. At the $X=5D$ position, the influence of tower drag diminishes, and the wake velocity begins to recover progressively, showing approximately 30% velocity increase at hub height. At the $X=7D$ position, compared with the $X=5D$ position, the wake velocity demonstrates an additional 20% increase in recovery magnitude.

Figure 5. Contours of vorticity.

3.2. Power Generations

Figure 7 shows the distribution of wind turbine power under different wind shear exponents. As the wind shear exponent increases from 0.12 to 0.15, there is no significant rise in power; instead, a slight decrease is observed. However, when the exponent increases from 0.15 to 0.22, the turbine’s power output sees a moderate increase. As the wind shear exponent further rises from 0.22 to 0.30, the rate of power growth becomes noticeably faster.

Figure 6. Vertical profiles of dimensionless axial velocity at different downstream positions under various wind shear indices.

Figure 7. Comparison of power at different wind shear exponents.

4. Conclusions

This study, using ALM and LES methods, explores the effects of yaw angles on the wake characteristics and power output of tandem wind turbines. The key conclusions are as follows:

Higher wind shear exponents result in faster wake recovery and increased spacing between tip vortices.

The expansion of tip vortices weakens as the wind shear exponent increases, with vortex breakdown occurring around two rotor diameters downstream. The tip vortices show an “elongated top, shortened bottom” pattern.

When the wind shear exponent increases from 0.12 to 0.15, wind turbine power experiences a slight decrease. However, as the exponent continues to rise, the overall energy extracted from the incoming flow increases. At a wind shear exponent of 0.30, power output is approximately 3% higher than at 0.12, corresponding to an absolute power increase of 613 KW.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Qian, Y., Wang, T., Yuan, Y. and Zhang, Y. (2020) Comparative Study on Wind Turbine Wakes Using a Modified Partially-Averaged Navier-Stokes Method and Large Eddy Simulation. Energy, 206, Article 118147.[CrossRef]
[2]	Wang, Z., Qiao, D., Tang, G., Wang, B., Yan, J. and Ou, J. (2022) An Identification Method of Floating Wind Turbine Tower Responses Using Deep Learning Technology in the Monitoring System. Ocean Engineering, 261, Article 112105.[CrossRef]
[3]	Fang, Y., Li, G., Duan, L., Han, Z. and Zhao, Y. (2021) Effect of Surge Motion on Rotor Aerodynamics and Wake Characteristics of a Floating Horizontal-Axis Wind Turbine. Energy, 218, Article 119519.[CrossRef]
[4]	Martínez-Tossas, L.A., Churchfield, M.J. and Leonardi, S. (2014) Large Eddy Simulations of the Flow Past Wind Turbines: Actuator Line and Disk Modeling. Wind Energy, 18, 1047-1060.[CrossRef]
[5]	Kim, T., Oh, S. and Yee, K. (2015) Improved Actuator Surface Method for Wind Turbine Application. Renewable Energy, 76, 16-26.[CrossRef]
[6]	Shen, W.Z., Zhu, W.J. and Sørensen, J.N. (2011) Actuator Line/Navier-Stokes Computations for the MEXICO Rotor: Comparison with Detailed Measurements. Wind Energy, 15, 811-825.[CrossRef]
[7]	Zhang, M., Arendshorst, M.G. and Stevens, R.J.A.M. (2018) Large Eddy Simulations of the Effect of Vertical Staggering in Large Wind Farms. Wind Energy, 22, 189-204.[CrossRef]
[8]	Onel, H.C. and Tuncer, I.H. (2021) Investigation of Wind Turbine Wakes and Wake Recovery in a Tandem Configuration Using Actuator Line Model with Les. Computers & Fluids, 220, Article 104872.[CrossRef]
[9]	Tu, Y., Zhang, K., Han, Z., Zhou, D. and Bilgen, O. (2023) Aerodynamic Characterization of Two Tandem Wind Turbines under Yaw Misalignment Control Using Actuator Line Model. Ocean Engineering, 281, Article 114992.[CrossRef]
[10]	Parada, L., Herrera, C., Flores, P. and Parada, V. (2018) Assessing the Energy Benefit of Using a Wind Turbine Micro-Siting Model. Renewable Energy, 118, 591-601.[CrossRef]
[11]	Zhu, X., Sun, C., Ouyang, H. and Du, Z. (2022) Numerical Investigation of the Effect of Towers and Nacelles on the Near Wake of a Horizontal-Axis Wind Turbine Model. Energy, 238, Article 121782.[CrossRef]