A Nonlinear Wind Turbine Wake Expansion Model Considering Atmospheric Stability and Ground Effects

Abstract

This study investigates the influence of atmospheric stability and ground effects on wind turbine wake recovery, challenging the conventional linear relationship between turbulence intensity and wake expansion coefficient. Through comprehensive field measurements and numerical simulations, we demonstrate that the linear wake expansion assumption is invalid at far-wake locations under high turbulence conditions, primarily due to ground effects. We propose a novel nonlinear wake expansion model that incorporates both atmospheric stability and ground effects by introducing a logarithmic relationship between the wake expansion coefficient and turbulence intensity. Validation results reveal the superior prediction accuracy of the proposed model compared to typical engineering wake models, with root mean square errors of wake wind speed predictions ranging from 0.04 to 0.063. The proposed model offers significant potential for optimizing wind farm layouts and enhancing overall wind energy production efficiency.

1. Introduction

Wind turbine wakes, characterized by reduced velocity and increased turbulence, significantly affect downstream turbine performance. These wakes are influenced by atmospheric conditions, ground effects, and turbine characteristics. Engineering wake models, which are crucial for wind farm design and operation, have evolved from simple linear models to complex Gaussian models. However, most existing models assume a linear relationship between wake expansion and turbulence intensity, often neglecting the effects of atmospheric stability and ground interactions. This assumption may be invalid under high turbulence conditions or in far-wake regions, especially in unstable atmospheric states where ground effects become significant. Therefore, there is a need for more sophisticated models that can accurately capture these complex wake behaviors.

Numerous wind tunnel [1,2,3] and field experiments [4,5,6] have demonstrated that atmospheric stability can significantly affect wind turbine wakes. Chamorro et al. [1] conducted an analysis of wind tunnel experiment data and found that under stable conditions, larger vertical gradients in wind speed compared to neutral conditions could enhance turbulence intensity in the wake region and extend its influence. Zhang et al. [3] performed wind tunnel experiments and observed that under unstable conditions, higher turbulence intensity compared to neutral conditions promotes wake recovery, resulting in a reduction in wake losses of approximately 15% and an increase in maximum turbulence intensity by 20%. Currently, engineering wake models that consider the influence of atmospheric stability primarily include the improved Jensen model by Peña et al. [7] and the wake model proposed by Cheng et al. [8]. The improved Jensen model predicts a wake profile resembling a “top hat”, which differs significantly from the measured Gaussian-like distribution shape. The Cheng model [8] has only been validated using numerical simulation results and lacks validation with field measurements, and thus its reliability requires further investigation.

Ground effects on wind turbine wakes are multifaceted and significant. The ground surface plays a crucial role in several aspects of wake behavior. These include the following: (1) Ground-induced turbulence and anisotropy: the ground contributes significantly to turbulence in the lower atmosphere. Under neutral conditions, a larger ground roughness increases the turbulence intensity, thereby accelerating wake recovery [9]. This ground-induced turbulence is also anisotropic, causing disparities between vertical and horizontal wake recovery rates [10]. (2) Wake profile shearing: the presence of the ground causes shearing of the wake profiles, introducing additional complexity to wakes. Several three-dimensional wake models have been proposed in recent years [11,12,13]. (3) Momentum interaction: in far-wake regions, the ground interaction of the wake impedes wake expansion toward the surface and may break down the linear wake expansion assumption [14]. At near and intermediate wake locations (approximately 3–5 rotor diameters downstream), the vertical distribution of the streamwise velocity deficit exhibits symmetry around the hub height. However, this symmetry breaks down in the far-wake region [15]. For simplicity, we focus on examining how ground effects influence wake expansion in this study.

This study addresses these issues by using numerical simulations to investigate the inhibitory effect of the ground on the wake expansion coefficient of wind turbines. A nonlinear wake expansion model that incorporates atmospheric stability and ground effects is proposed based on both measured data and relevant numerical simulation results. The remainder of this paper is organized as follows. The effects of ground and atmospheric stability on wind turbine wake expansion are studied in Section 2. The improved Jensen wake model and several typical models are introduced in Section 3. A nonlinear wake expansion model based on the experimental and numerical results is proposed in Section 4 and validated in Section 5. The conclusions are presented in Section 6.

2. Effects of the Ground on the Wind Turbine Wake Expansion

2.1. FullRF Turbulence Model for Wake Modeling

The study applies the FullRF turbulence model [16] and the actuator disk model to simulate the wind turbine wake, with the control equations adopted as follows:

$${\displaystyle \frac{\partial }{\partial x_i}}\left(\rho U_i\right)=0$$

(1)

$${\displaystyle \frac{\partial }{\partial t}}(\rho U_i)+{\displaystyle \frac{\partial }{\partial x_j}}(\rho U_i{U}_j)=-{\displaystyle \frac{\partial p}{\partial x_i}}+{\displaystyle \frac{\partial }{\partial x_j}}\left[(\mu +{\mu }_t)\left({\displaystyle \frac{\partial U_i}{\partial x_j}}+{\displaystyle \frac{\partial U_j}{\partial x_i}}\right)\right]+S_{u,i}$$

(2)

$${\displaystyle \frac{\partial }{\partial t}}(\rho \Theta )+{\displaystyle \frac{\partial }{\partial x_i}}(\rho U_i\Theta )={\displaystyle \frac{\partial }{\partial x_i}}\left[\left({\displaystyle \frac{\mu }{Pr}}+{\alpha }_t\right){\displaystyle \frac{\partial \Theta }{\partial x_i}}\right]$$

(3)

$${\displaystyle \frac{\partial }{\partial t}}(\rho k)+{\displaystyle \frac{\partial }{\partial x_i}}(\rho U_{i}k)={\displaystyle \frac{\partial }{\partial x_j}}\left[\left(\mu +{\displaystyle \frac{{\mu }_t}{{\sigma }_k}}\right){\displaystyle \frac{\partial k}{\partial x_j}}\right]+\mathcal{P}+\mathcal{B}-\rho \epsilon -S_{k,\mathrm{ASL}}$$

(4)

$${\displaystyle \frac{\partial }{\partial t}}(\rho \epsilon )+{\displaystyle \frac{\partial }{\partial x_i}}(\rho U_i\epsilon )={\displaystyle \frac{\partial }{\partial x_j}}\left[\left(\mu +{\displaystyle \frac{{\mu }_t}{{\sigma }_{\epsilon }}}\right){\displaystyle \frac{\partial \epsilon }{\partial x_j}}\right]+\left(C_{\epsilon 1}\mathcal{P}-\rho C_{\epsilon 2}\epsilon +C_{\epsilon 3}\mathcal{B}\right){\displaystyle \frac{\epsilon }{k}}+S_{\epsilon ,\mathrm{wake}}$$

(5)

where $U_i$ and $U_j$ represent the velocity components along the $x_i$ and $x_j$ axes, respectively, p is pressure, $\mu$ is molecular viscosity, µ_t is the turbulent viscosity, and $S_{u,i}$ is the momentum source term exerted by the wind turbine on the $x_i$ axis. $\Theta$ is the potential temperature, $Pr=0.9$ is the laminar Prandtl number, α_t ≡ µ_t/Pr_t is the turbulent heat conductivity, and $P{r}_t\equiv {\varphi }_h/{\varphi }_m$ is the turbulent Prandtl number. The similarity functions ${\varphi }_m$ and ${\varphi }_h$ are described in Section 2.1.1. k is the turbulence kinetic energy (TKE), ε is the TKE dissipation, $C_{\mu }=0.033$, ${\sigma }_k=1.0$, ${\sigma }_{\epsilon }=1.3$, and $C_{\epsilon 2}=1.92$. $S_{k,\mathrm{ASL}}$ and $S_{\epsilon ,\mathrm{wake}}$ represent the source terms for turbulence kinetic energy and turbulence kinetic energy dissipation rate, and are also described in Section 2.1.1. The turbulent kinetic energy generation rate $\mathcal{P}$ due to the mean velocity gradient and the turbulent kinetic energy generation rate $\mathcal{B}$ due to the buoyancy force can be expressed under the Boussinesq assumption,

$$\mathcal{P}\equiv -\rho \overline{u_i^{\prime }{u}_j^{\prime }}{\displaystyle \frac{\partial U_i}{\partial x_j}}={\mu }_t\left({\displaystyle \frac{\partial U_i}{\partial x_j}}+{\displaystyle \frac{\partial U_j}{\partial x_i}}\right){\displaystyle \frac{\partial U_i}{\partial x_j}}$$

(6)

$$\mathcal{B}\equiv {\displaystyle \frac{g_j}{{\Theta }_0}}\overline{u_j^{\prime }{\theta }^{\prime }}=-{\displaystyle \frac{g_j}{{\Theta }_0}}{\displaystyle \frac{{\mu }_t}{P{r}_t}}{\displaystyle \frac{\partial \Theta }{\partial x_j}}$$

(7)

where $u_j^{\prime }$ and $g_j$ are the components of wind speed fluctuation and gravitational acceleration along the $x_j$ axis, and ${\theta }^{\prime }$ represents the potential temperature fluctuation.

This section uses the BEM method to simulate the wind turbine rotor. In the BEM-based actuator disk model, the rotor plane consists of N actuator lines, each of which is split into M-element sections (Figure 1). The element section collects the local velocity and rotor speed $\Omega$ to calculate the element force, and applies this force to the neighboring cells of the element section. The reference velocity is first assessed from the disk-averaged velocity and is then applied to evaluate the rotor speed.

When transforming the local velocity at the blade element into polar velocity components ($\Omega r$, $u_{\theta }$, and $u_n$), the force of the blade element is

$$\overrightarrow{\Delta F}=\rho {\displaystyle \frac{Bc\Delta \theta \Delta r}{4\pi }}\left(C_L\overrightarrow{e_L}+C_D\overrightarrow{e_D}\right)\left[u_n^2+{\left(\Omega r+u_{\theta }\right)}^2\right]$$

(8)

where $B$ is the number of blades, $c$ is the chord length, and ∆r is the element section length. The drag coefficient of the element section, $C_D$, and its lift coefficient, $C_L$, which are functions of the attack angle α, are estimated from XFOIL [17] and then corrected by the three-dimensional rotational effects of the blades based on Du et al. [18]. According to Figure 1, α = φ − (β + γ), where φ = arctan [u_n/(Ωr + u_θ)] is the flow angle, β is the blade installation angle, and γ is the pitch angle.

The element force is distributed across neighboring cells. The force added to a cell is calculated by

$$\overrightarrow{\Delta F_{\mathrm{cell}}}={\displaystyle \sum _{j=0}^{N-1}{\displaystyle \sum _{i=0}^{M-1}\frac{1}{s^3{\pi }^{3/2}}\exp \left(-\frac{s_{i,j}^2}{s^2}\right)}}\overrightarrow{\Delta F_{i,j}}{f}_i^{\mathrm{tip}}{f}_i^{\mathrm{hub}}\Delta V_{\mathrm{cell}}$$

(9)

where s_i,j is the distance of the i-th element to the cell and s is the cut-off length scale that takes a value between two and three cell sizes [19]. $f_i^{\mathrm{tip}}$ and $f_i^{\mathrm{hub}}$ are the Prandtl tip loss and hub loss functions [20]:

$$f_i^{\mathrm{tip}}={\displaystyle \frac{2}{\pi }}\arccos \left[\exp \left({\displaystyle \frac{B\left(R-r_i\right)}{2r_i\sin {\varphi }_i}}\right)\right]$$

(10)

$$f_i^{\mathrm{hub}}={\displaystyle \frac{2}{\pi }}\arccos \left[\exp \left({\displaystyle \frac{B\left(r_i-R_{\mathrm{hub}}\right)}{2r_i\sin {\varphi }_i}}\right)\right]$$

(11)

where R is the rotor radius, R_hub is the hub radius, and r is the radial distance between the element and the rotor center.

2.1.1. Turbulence Modeling

In the FullRF turbulence model, a TKE source term and the coefficient $C_{\epsilon 3}$ are calibrated to keep flow homogeneity:

$$S_{k,\mathrm{ASL}}=[{\varphi }_m-{\varphi }_{\epsilon }-\zeta +{\displaystyle \frac{C_{\mu }^{-1/2}{\kappa }^2}{{\sigma }_k}}\frac{{\zeta }^2}{{\varphi }_m}({\varphi }_k^{″}-\frac{{\varphi }_k^{\prime }{\varphi }_m^{\prime }}{{\varphi }_m}+\frac{{\varphi }_k^{\prime }}{\zeta }\left\}\right]\frac{u_{*}^3}{\kappa z}$$

(12)

$$C_{\epsilon 3}={\displaystyle \frac{C_{\epsilon 1}{\varphi }_m-C_{\epsilon 2}{\varphi }_{\epsilon }}{\zeta }}+{\displaystyle \frac{C_{\mu }^{-1/2}{\kappa }^2}{{\sigma }_{\epsilon }}}{\displaystyle \frac{{\varphi }_k}{{\varphi }_m}}(\frac{\zeta {\varphi }_{\epsilon }^{″}}{{\varphi }_{\epsilon }}-\frac{\zeta {\varphi }_{\epsilon }^{\prime }{\varphi }_m^{\prime }}{{\varphi }_{\epsilon }{\varphi }_m}-\frac{{\varphi }_{\epsilon }^{\prime }}{{\varphi }_{\epsilon }}+\frac{{\varphi }_m^{\prime }}{{\varphi }_m}+{\displaystyle \frac{1}{\zeta }})$$

(13)

where the atmospheric stability parameter $\zeta =z/L$, z is the height above the ground, and $L$ is the Monin–Obukhov length:

$$L\equiv -\frac{u_{*}}{\kappa {\displaystyle \frac{g}{\Phi }}\overline{w^{\prime }{\theta }^{\prime }}}$$

(14)

in which $\kappa$ is the von Kármán constant (0.4), g is the acceleration due to gravity (9.8 m/s²), and $u_{*}=\sqrt{-\overline{u^{\prime }{w}^{\prime }}}$ is the friction velocity [16].

The similarity functions for the FullRF model are

$${\varphi }_m(\zeta )=\left\{\begin{array}{ll}{(1-16\zeta )}^{-1/4}\hfill & \zeta <0\hfill \\ {(1+40\zeta )}^{1/4}\hfill & 0<\zeta \hfill \end{array}\right.$$

(15)

$${\varphi }_h(\zeta )=\left\{\begin{array}{ll}0.9{(1-16\zeta )}^{-1/2}\hfill & -2<\zeta <0\hfill \\ 0.9{\left(1+5\zeta \right)}^{1/4}\hfill & 0<\zeta <1\hfill \end{array}\right.$$

(16)

$${\varphi }_{\epsilon }(\zeta )\equiv {\displaystyle \frac{\kappa z}{u_{*}^3}}\epsilon =\left\{\begin{array}{ll}1-\zeta \hfill & \zeta <0\hfill \\ {\varphi }_m-\zeta \hfill & \zeta >0\hfill \end{array}\right.$$

(17)

$${\varphi }_k(\zeta )\equiv {\displaystyle \frac{\sqrt{C_{\mu }}}{u_{*}^2}}k=\sqrt{{\displaystyle \frac{{\varphi }_{\epsilon }(\zeta )}{{\varphi }_m(\zeta )}}}$$

(18)

The source term $S_{\epsilon ,\mathrm{wake}}$ for turbulent kinetic energy dissipation is applied within a cylindrical region downstream of the wind turbine with a distance of 0.25 times the rotor diameter, aiming to correct the issue of rapid wake recovery induced by the $k-\epsilon$ standard turbulence model.

$$S_{\epsilon ,\mathrm{wake}}=\rho C_{\epsilon 4}{\displaystyle \frac{{\mathcal{P}}^2}{\rho k}}$$

(19)

in which $C_{\epsilon 4}=0.37$.

In addition, Alinot and Masson [21] and M.P. van der Laan [22] proposed two near-surface turbulence models based on the Businger–Dyer similarity function, which are referred to as the AM model and the Laan model. Section 2.2.3 will study the effects of atmospheric stability on wind turbine wake for the aforementioned models.

2.1.2. Boundary Conditions

The boundary conditions consistent with similarity functions are applied to model the atmospheric boundary stratification. The inlet profiles of wind speed, potential temperature, TKE, and its dissipation are given by

$$U\left(z\right)={\displaystyle {\int }_{z_0}^z\frac{u_{*}}{\kappa z}{\varphi }_m\left(\frac{z}{L}\right)\mathrm{d}z}$$

(20)

$$\Theta \left(z\right)={\Theta }_0+{\displaystyle {\int }_{z_0}^z\frac{{\theta }_{*}}{\kappa z}{\varphi }_h\left(\frac{z}{L}\right)\mathrm{d}z}$$

(21)

$$\epsilon \left(z\right)={\displaystyle \frac{u_{*}^3}{\kappa z}}{\varphi }_{\epsilon }\left({\displaystyle \frac{z}{L}}\right)$$

(22)

$$k\left(z\right)={\displaystyle \frac{u_{*}^2}{\sqrt{C_{\mu }}}}{\varphi }_k\left({\displaystyle \frac{z}{L}}\right)$$

(23)

where z₀ is the aerodynamic roughness length and ${\theta }_{*}=-\overline{w^{\prime }{\theta }^{\prime }}/u_{*}$ is the scaling temperature.

The vertical profiles of wind speed are estimated in numerical integration. Zero gradients of $U$, $\Theta$, $\epsilon$, and $k$ are applied at the outlet. For the top boundary, the upstream flow properties are maintained at constant. The left and right sides of the computational domain are symmetrical. The near-wall treatment of Temel et al. [23] is implemented in the near-ground region to calculate the turbulent dissipation rate and TKE production $G_{k,p}$ owing to shear and buoyancy:

$${\epsilon }_p={\displaystyle \frac{u_{*k}^3}{\kappa z_p}}{\varphi }_{\epsilon ,p}$$

(24)

$$G_{k,p}={\displaystyle \frac{{\tau }_w^2}{\rho \kappa u_{*k}\left(z_p+z_0\right)}}-{\displaystyle \frac{{\theta }_{*}\left|g\right|}{{\Theta }_0\left(z_p+z_0\right)}}$$

(25)

where the subscript p denotes the first cell center above the ground, the equivalent friction wind speed is $u_{*k}=C_{\mu }^{0.25}{\varphi }_{k,p}^{-0.5}{k}_p^{0.5}$, and the wall shear stress is ${\tau }_w={\mu }_{t,p}\left(\mathrm{d}U/\mathrm{d}z\right)$. The eddy viscosity ${\mu }_{t,p}$ and the turbulent heat conductivity are imposed at the first cell center above the ground, as in [24]:

$${\mu }_{t,p}={\displaystyle \frac{u_{*k}{z}_p}{{\displaystyle {\int }_{z_0}^{z_p}\frac{1}{\kappa z}{\varphi }_m\left(\frac{z}{L}\right)\mathrm{d}z}}}-\mu$$

(26)

$${\alpha }_{t,p}={\mu }_{t,p}{\varphi }_{m,p}/{\varphi }_{h,p}-\mu /Pr$$

(27)

2.2. Model Validation

2.2.1. Test Case

To validate the FullRF model, the scenarios outlined in reference [16] only involved wind turbines with capacities below 1 MW and lacked examination under neutral and unstable atmospheric conditions. In this study, the wake measurement data from the Jingbian wind farm [4] are used to further substantiate the FullRF model. Figure 2 shows the wind farm’s topography and the layout of wind turbines; the X and Y axes are oriented towards the east and north, respectively. The southern region of the wind farm, characterized by a valley with complex topography, contrasts with the relatively flat northern expanse. Experimental wind turbine #14 has a rated capacity of 2 MW, rotor diameter of 90 m, and hub height of 67 m, and is situated in a transitional zone between the valley and flat terrain. Consequently, wake measurements were executed using masts M1 and M3, positioned to the south and north of turbine #14, at horizontal distances approximately 1.45 times and 2.15 times the rotor diameter, respectively. The masts were equipped with Thies First Class cup anemometers, temperature and wind direction sensors, and Metek 3D sonic anemometers, which provided wind speed and temperature data at a frequency of 35 Hz. This study collected approximately 190 days of sonic anemometer measurements and 310 days of measurement data from masts and wind farm SCADA systems. Since the atmospheric boundary layer classification relies on data from the sonic anemometers, this study uses 190 days of data for analysis.

Masts were strategically positioned on both the north and south sides of wind turbine #14, as shown in Figure 3. On each mast, anemometers were installed at heights of 30 m, 50 m, 60 m, and 70 m above the ground. This setup enabled the acquisition of the wake wind velocity and turbulence intensity profiles at axial distances of 1.45D and 2.15D behind the turbine. Given that wind turbines are typically spaced more than 4D apart, calculations of wind resources and forecasts of wind power are primarily concerned with the far-wake region of turbines. In addition, the nacelle wind speed, when combined with the nacelle transfer function, provides a partial estimation of the wind speed at the turbine. To augment these data, an experimental setup involving an additional turbine, #12, was introduced near the M1 mast to capture the wake profile at a further axial distance of 5D, as shown in Figure 4.

The wind speed of the hub corresponding to the wake measurement data of the wind turbine is 6 m/s, the aerodynamic roughness length is 0.05 m, and the Monin–Obukhov length $L$ is estimated in every 10 min according to Equation (14) based on the sonic anemometer measurements at 30 m above the ground [16]. The flow parameters corresponding to each stability regime are listed in

Table 1. Inflow parameters for wake simulations of the Haizhuang 2 WM wind turbine at the hub height.

Atmosphere Stability	$\mathit{L}$ (m)	$\mathbf{Measured} {\mathit{I}}_{\mathit{u}}$ (%)	$\mathbf{CFD} {\mathit{u}}_{*}$ (m/s)		$\mathbf{CFD} {\mathit{I}}_{\mathit{u}}$ (%)
			Classical	FullRF	Classical	FullRF
Unstable	−30	17.6	0.425		22.8
Neutral	$\infty$	11.6	0.333		10.6
Stable	30	7.8	0.131	0.208	4	9

, where the turbulence intensity $I_u$ is measured or estimated at the hub height of 67 m and “Classical” stands for the Businger–Dyer similarity functions [25,26].

2.2.2. Computational Domain, Meshing, and Solver Settings

The calculation domains of the wind turbine wake simulation are 20D, 10D, and 10D, respectively (Figure 5). To effectively incorporate the wind turbine source term into the computational cells and capture the wake structure, the mesh of the computational domain is segmented into four refinement levels. A cell at level i can be subdivided into four cells at level i + 1. The background grid with refinement level 0 contains 100 (long) × 60 (wide) × 60 (high) cells. The background grid is uniformly divided in the horizontal plane, and refined near the ground in the height direction, and the height of the grid is set to 7.38z₀ [27]. The grid refinement level in the region of the wind turbine actuator is set to 3, distributing approximately 80 grid nodes across the diameter and length of the turbine [28]. The actuator disk regions 5D and 10D downstream are refined to levels 2 and 3, respectively, to ensure detailed capture of the wake structures. The computational domain contains approximately 1.6 million cells. In this study, turbulence models are implemented in OpenFOAM [29] using the finite-volume method. A large time-step transient solver using the PIMPLE algorithm is applied to simulate wakes of wind turbines. During the iteration, the source terms $S_{u,i}$ and $S_{\epsilon ,\mathrm{wake}}$ are modeled based on local flow information via user-specified finite volume options. For temporal and spatial discretization in the simulations, a second-order backward difference scheme and a second-order central difference scheme are applied, respectively.

2.2.3. Results

Figure 6 presents a comparison between the relative wind speeds predicted by the wake model at various axial distances and the corresponding empirical measurements. Here, FullRF denotes the utilization of the FullRF turbulence model coupled with the BEM-based actuator disk model, whereas FullRF-CT refers to the adoption of the FullRF turbulence model integrated with an actuator disk model that employs thrust coefficients.

Figure 7 illustrates the correlations between the root-mean-square error (RMSE) of wind speed at different axial distances. Since the actuator disk model based on the thrust coefficient does not add a momentum source term to the nacelle position, the actuating disk resembles a ring shape, forming a bimodal phenomenon when the wind speed decreases at a hub height of 1.45D after the wind wheel. Compared with the AM model and the Laan model, the wake predicted by the FullRF and FullRF-CT models has a lower RMSE relative to the wind speed, and is in good agreement with the measured values. Under stable conditions, the AM model and the Laan models significantly overestimate the wind speed loss in the wake region, and overestimate the wake loss by approximately 40% at 5D, and the corresponding RMSE is maintained at approximately 0.15. The RMSE of the wake predicted by the FullRF and FullRF-CT models is only about 0.7 on average, which effectively improves the phenomenon of overestimating the wind speed loss caused by the Businger–Dyer similarity functions. Under unstable conditions, the turbulence intensity of the reference flow provided by the model is 22.8%, which is greater than the measured value of 18%. This accelerates the recovery speed of the wake in the simulation results, resulting in the wind speed predicted by the model being higher than the measured values (Figure 6).

2.3. Effects of the Ground on Wake Expansion

Figure 8 depicts the distribution of wake wind speeds in the vertical cross-section of a wind turbine at various distances downstream under unstable, neutral, and stable atmospheric conditions. The dark line encircling the wake indicates the wake boundary, which is characterized by a wake deficit of 0.05. The illustration reveals varying degrees of downstream wake center displacement at the 10D position behind the turbine across different operating conditions. The displacement is most pronounced under unstable conditions, which feature high turbulence intensity. Under such conditions, wake expansion occurs more rapidly compared to scenarios with lower turbulence. Additionally, the interaction between the ground and the wake is enhanced under high turbulence. This will cause the wake’s centerline to shift downward towards the ground in the far-wake region, and the linear expansion of the wake at hub height will be suppressed. Further investigations are conducted in Section 4 to assess whether the wake maintains a linear expansion under the influence of high turbulence and ground effects.

3. Typical Engineering Wake Models with Linear Wake Expansion

3.1. A Modified Jensen Model Considering Atmospheric Stability Conditions

In the Jensen wake model [30], the wind speed deficit at a distance x downstream of the rotor is calculated as follows:

$${\displaystyle \frac{\Delta U(x,r)}{U_{\infty }}}={\displaystyle \frac{1-\sqrt{1-C_T}}{{(1+2k_{w}x/D)}^2}}$$

(28)

where $U_{\infty }$ is the inflow velocity, $C_T$ is the thrust coefficient, and $k_w$ is the wake expansion rate and could be used to determine the wake boundary; at $r>D/2+k_{w}x$, $\Delta U(x)=0$ where $r=y-y_c$ represents the radial distance. In this paper, the coordinates of the wind turbine rotor’s rotation center are set as the origin, i.e., $x_c=y_c=z_c=0$.

For wind turbine wakes in a neutral atmospheric boundary layer, the suggested values of k_w are 0.075 for onshore cases and 0.04 or 0.05 for offshore ones [31,32]. To account for the influence of atmospheric stability on the wake, Peña et al. [33] related the wake expansion coefficient to the hub-height incoming wind speed U_H and the friction velocity $u_{*}$:

$$k_w=u_{*}/U_H$$

(29)

$$u_{*}={\displaystyle \frac{U_H}{{\displaystyle {\int }_{z_0}^H\frac{1}{z}}{\varphi }_m(\zeta )\mathrm{d}z}}$$

(30)

where H represents the hub height.

3.2. Guassian-Shaped Wake Models

Because of the assumption of a top hat distribution for wake deficit profiles, the modified Jensen model cannot capture the radial dependence of the wake. In fact, the wake deficit has an approximately Gaussian symmetrical shape after some downwind distances [1]. By applying the conservation of mass and momentum, Bastankhah and Porté-Agel [34] suggested replacing the top hat assumption with a Gaussian distribution for the wake deficit in the wake:

$${\displaystyle \frac{\Delta U(x,r)}{U_{\infty }}}=\left[1-\sqrt{1-{\displaystyle \frac{C_T}{8{(\sigma /D)}^2}}}\right]\exp \left(-{\displaystyle \frac{r^2}{2{\sigma }^2}}\right)$$

(31)

where $\sigma$ is the wake width.

As shown in various numerical and experimental studies on wind turbine wakes [34,35], the wake width is approximately linear with x after some downstream distance:

$${\displaystyle \frac{\sigma }{D}}=k^{*}{\displaystyle \frac{x-x_0}{D}}+{\sigma }_0\equiv k^{*}{\displaystyle \frac{x}{D}}+{\sigma }_R$$

(32)

in which $k^{*}\equiv \partial \sigma /\partial x$ represents the wake expansion coefficient corresponding to the Gaussian distribution wake model and $x_0$ characterizes the near-wake length, defined as the axial distance at which the wind speed deficit exhibits a Gaussian distribution-like profile. ${\sigma }_R$ is the parameter equivalent to the wake width at (x = 0), as assumed in the Gaussian distribution wake model. Since the wind speed deficit at the rotor varies from the Gaussian distribution, $k^{*}$ is often considered a linear function of the longitudinal turbulence intensity [34]. ${\sigma }_R$ is then divided into two categories, one being a linear function of $k^{*}$ [34] and the other related to $k^{*}$, as well as the near-wake length x₀ [36]:

$${\sigma }_R=-{\displaystyle \frac{x_0}{D}}{k}^{*}+{\sigma }_0$$

(33)

Bastankhah et al. [36] divided the wake of a wind turbine into three regions: the near-wake core region, the outer atmospheric free-flow region, and the boundary layer region between them (Figure 9). They derived a formula for calculation by relating the growth rate of the boundary layer thickness to the turbulence intensity and the difference in wind speeds inside and outside the wake. The formula used is as follows:

$${\displaystyle \frac{x_0}{D}}={\displaystyle \frac{1+\sqrt{1-C_T}}{\sqrt{2}[{\alpha }_{*}{I}_n+{\beta }_{*}(1-\sqrt{{1-C}_T})]}}$$

(34)

with ${\beta }_{*}=0.154$, and ${\alpha }_{*}$ is typically set to 2.32 in wind tunnel experiments. Based on the analysis of field wake observations by Fuertes et al. [37], the near-wake length $x_0$ of wind turbines in wind farms is observed to be smaller than the values observed in wind tunnels under equivalent turbulence intensity and thrust coefficient conditions. Therefore, in wind farms, ${\alpha }_{*}$ is typically set to the value recommended by Fuertes et al. [37], which is 3.6. ${\sigma }_0=1/\sqrt{8}$ for unyawed conditions.

Turbulence plays a significant role in the evolution of wind turbine wakes and is generally considered to have a linear relationship with the wake expansion coefficient $k^{*}$. This study introduces three typical Gaussian wake models: the Fuertes model [37], Cheng model [8], and Campagnolo model [38] (

Table 2. Typical Gaussian wake models, where $I_u$ and $I_v$ are longitudinal turbulence intensity and lateral turbulence intensities, respectively.

Model	Fuertes	Cheng	Campagnolo
Scale	The experimental measurement of the nacelle lidar	Considering lateral turbulence	Wind tunnel experiment
$k^{*}$	$0.35I_u$	$0.223I_v+0.022$	$0.089I_u+0.027$
${\sigma }_R$	$-1.91k^{*}+0.34$		Equation (33), ${\alpha }_{*}=0.952$, ${\beta }_{*}=0.262$

). These models are compared with the engineering model proposed in this study for predicting wind turbine wake wind speeds under different atmospheric stability conditions.

3.3. Estimation of the Streamwise Turbulence Intensity at Hub Height

To assess the prediction accuracy of the model, this study first utilizes existing measurements of longitudinal turbulence intensity. If such measurements are unavailable, a proposed similarity function is employed for estimation. Based on the similarity function ${\sigma }_u/u_{*}$ at a height of 70 m,

$$I_u\equiv {\displaystyle \frac{{\sigma }_u}{U_H}}={\displaystyle \frac{2.24{\left(1+b_0\zeta \right)}^{1/2}}{{\displaystyle {\int }_{z_0}^H}\kappa z{\varphi }_m\left(\zeta \right)\mathrm{d}z}}$$

(35)

the FullRF similarity function is employed for ${\varphi }_m$, with $b_0$ set to 0.4 and −0.5 for stable and unstable conditions, respectively.

The lateral turbulence intensity is estimated using the method recommended in the ESDU [39]:

$$I_v=\left[1-0.22{\cos }^4\left({\displaystyle \frac{\pi H}{2h_{\mathrm{ABL}}}}\right)\right]I_u$$

(36)

in which $h_{\mathrm{ABL}}=u_{*}/(12\Omega \sin \varphi )$ represents the boundary layer height, Ω = 72.9 × 10⁻⁶ rad/s denotes the Earth’s rotational angular velocity, and ϕ represents the local latitude.

4. The Proposed Engineering Wake Expansion Model

Atmospheric stability significantly influences the dynamics of wind turbine wakes by affecting turbulence intensity. Consequently, this study delineates the impact of atmospheric stability by establishing a correlation between turbulence intensity and the wake expansion coefficient. Equation (33) provides a method for estimating the wake expansion width based on the length of the near wake. This study predominantly utilizes numerical simulations to construct the model under high turbulence scenarios. Furthermore, under stable conditions, classical similarity functions often overpredict the vertical gradients of wind speed and underpredict the incoming turbulence intensity [16]. To address these discrepancies, the chapter uses similarity functions, adjusted based on measured data from a wind farm in China, to simulate wake effects. The simulation results are then used to create a data training set and a validation set, which are employed to develop an engineering wake model that considers atmospheric stability. For specific details, refer to Figure 10.

4.1. Wake Model Development Data Set

The dataset consists of two parts: training data used to fit the correlation between $k^{*}$ and $I_u$, and testing data used to validate the reliability of the developed wake model. The dataset description is shown in

Table 3. Basic information of the wake dataset.

Experiment Type	Cases	Subcases	D (m)	H (m)	Wake Range (D)
Wind tunnel tests	Dou2019 [40,41]	TSR = 4, 5, 6	0.2	0.75	4.5~10
	WiRE-01 [42]	/	0.15	0.125	4~10
	Ruland-913 [43]	In turbulent flows	0.9	1.12	2.5~8.5
	G1 [38]	Three offshore cases	1.1	0.83	5~10
	Hancock2014 [2,44,45]	Neutral, unstable, and stable cases	0.416	0.3	3~10
Field experiments	Nibe-B [46,47]	$C_T$ = 0.67, 0.77, 0.82	40	45	2.5~7.5
	Liberty C96 [37]	/	96	80	0.6~10
	Vestas V80-2MW [9]	Neutral (LES): 4 types of z0 and three stability classes (LES)	80	70	3~15
	Danwin [48,49,50]	Neutral: $C_T$ = 0.65, 0.82 Non-neutral: (experiments + RANS)	23	35	4.2~9.6
	Nordtank [5]	Three stability classes (experiments + RANS)	41	36	2~5
	Haizhuang [4]	Three stability classes (experiments + RANS)	93	67	1.45, 2.15, 5
	Nibe-B [46,47]	$C_T$ = 0.67, 0.77,0.82	40	45	2.5~7.5

. The dataset includes wind speed data from the wake region of 13 types of wind turbines, obtained through wind tunnel experiments, field measurements, and numerical simulations, along with relevant inflow information. The shaded sections in the table indicate cases that involve non-neutral operating conditions. The cases are primarily matched with the name of the wind turbine, and when the turbine name is unknown, they correspond to the author and publication date of the data source. The subclass denotes the presence of multiple operating environments or conditions for the same wind turbine, where TSR represents the Tip Speed Ratio. The three stability levels indicate the presence of neutral, unstable, and stable operating environments. In the subclass, LES stands for Large Eddy Simulation, with inflow conditions set based on the Businger–Dyer similarity functions, considering the influence of the atmospheric boundary layer thickness. The effect of the atmospheric boundary layer thickness is confirmed by applying a pressure gradient within the computational domain to ensure that the Reynolds stress decreases to zero at the top of the atmospheric layer. RANS refers to Reynolds-Averaged Navier–Stokes simulations using the FullRF turbulence model to simulate wind turbine wakes.

The wind turbine operating conditions corresponding to wind tunnel experiments are as follows: (1) A two-bladed wind turbine in the Dou 2019 case [40,41], operating at different tip-speed ratios (TSR = 4, 5, 6), pitch angles, and yaw angles. (2) A model wind turbine (WiRE-01) with a diameter of 15 cm, operating at the optimum tip-speed ratio (${\lambda }_{opt}\approx 3.8$), corresponding to a thrust coefficient $C_T\approx 0.78$ [42]. (3) A six-bladed Rutland 913 model wind turbine operating under turbulent conditions ($I_u=14.5\%$) [43]. (4) The G1 turbine operates in two environments, one simulating offshore wind farm conditions with a turbulence intensity of 6.1% and the other simulating onshore wind farm conditions with a turbulence intensity of 11% [38]. (5) The Hancock 2014 case involves a three-bladed model wind turbine operating under neutral, stable, and unstable conditions [2,44,45].

For non-wind-tunnel experiments, the Vestas V80-2MW wind turbine only has wake results from LES results [9], while other turbines have field-measured data. Turbines with observed wake data and operating conditions include the following: (1) A 40 m diameter three-bladed Nibe-B 630 kW wind turbine operating at 33 rpm [46,47]. (2) The 2.5 MW Liberty C96 wind turbine with a thrust coefficient of approximately 0.82, capable of measuring wake wind speeds from 0.6 to 10 rotor diameters using a laser radar mounted on the nacelle [37]. (3) The three-bladed Danwin 180 kW wind turbine operating under different atmospheric thermal stability conditions [48,49,50]. (4) The actively stalled 500 kW Nordtank wind turbine equipped with a pulsed laser radar mounted on the nacelle, capable of measuring wake wind speeds behind the rotor at different atmospheric stability conditions [5]. (5) The Haizhuang 2 MW wind turbine, with wind measurement towers installed on both the north and south sides and related wind turbine SCADA operational data, providing wind speed data at positions 1.45D, 2.15D, and 5D behind the rotor under various stability conditions [4].

Table 4. Supplementary wind turbine wake simulation cases.

Turbine	${\mathit{z}}_0$ (m)	$\mathit{L}$ (m)	${\mathit{u}}_{*}$ (m)	${\mathit{I}}_{\mathit{u}}$ (%)	${\mathit{k}}^{*}$ (10−3)
Haizhuang	5 × 10−8	∞	0.114	3.6	5.79
	5 × 10−6		0.146	4.7	6.28
	5 × 10−5		0.17	5.4	8.07
	0.05	−100	0.382	16.2	37.47
		−1000	0.5122	17.4	44.05
	0.5	−100	0.601	25.6	69.38
		−50	0.655	31.6	70.13
		−20	0.755	45	83.43
Nordtank	0.2	∞	0.542	14.6	51.26
Danwin	5 × 10−4	35	0.198	4.5	3.23

provides supplementary parameters for the wake simulation cases of wind turbines, with the wake expansion coefficient obtained through the least squares fitting of numerical simulation data.

4.1.1. The Training Dataset

Table 5. The hub-height turbulence intensity and thrust coefficient of wind turbines.

Training Cases		${\mathit{I}}_{\mathit{u}}$ (%)	${\mathit{C}}_{\mathit{T}}$
Wind tunnel experiments	Dou2019 (TSR = 4)	1	0.85
	G1: offshore	6.1	0.79, 0.73, 0.68
	Rutland-913	14.5	0.94
Field measurements	Nibe-B	11	0.77, 0.82
	Liberty C96	1.6~17	0.82
LES simulations	Vestas: neutral	4.8~13.4	0.8
	Vestas: three stability classes	6.5~10	0.8
RANS simulations	Haizhuang	3.6~44.7	0.84
	Nordtank	6.1~18.3	0.83
	Danwin	4.5~10	0.82

presents turbulence intensity and thrust coefficients for the training data. The turbulence intensity provided in the training data ranges from 1.6% to 45%. Given that LES numerical results from Vestas wind turbines frequently informed the development of wake models, these simulations have also been incorporated into the training dataset for this study [51,52,53].

4.1.2. Validation Data

Apart from the training data, other cases will be employed to verify the accuracy and reliability of the proposed model in predicting wind turbine wakes. The inflow conditions and thrust coefficients for these validation cases are listed in

Table 6. The inflow information and thrust coefficient of wind turbines for the validation cases.

Case Validation		${\mathit{U}}_{\mathit{\infty }}$	$\mathit{L}$ (m)	${\mathit{I}}_{\mathit{u}}$ (%)	${\mathit{C}}_{\mathit{T}}$
Wind tunnel	Dou2019 (TSR = 5, 6)	6	∞	1	0.91, 0.94
	WIRE-01	5	∞	7
	Hancock2014	2.3, 2.3, 1.47	0.956, ∞, −1.26	8.5, 6.6, 5.3	0.42, 0.48, 0.48
Field	Nibe-B	11.52	∞	10.5	0.67
	Danwin	8, 11, 8	−50, ∞, 90.6	9.7, 6, 7.6	0.82, 0.65, 0.82
	Nordtank	6.82, 7.03, 6.76	−84.8, ∞, 29	14, 15, 10	0.71,0.75,0.83

.

4.2. The LogIu Engineering Wake Model

The relationship between the wake expansion coefficient $k^{*}$ and the longitudinal turbulence intensity $I_u$ obtained from the training data is shown in Figure 11. The training data suggest a proportional relationship between the wake expansion coefficient and the longitudinal turbulence when the longitudinal turbulence intensity is below 10%. Beyond this threshold, the rate of change of the wake expansion coefficient diminishes with increasing turbulence intensity. In light of supplementary RANS numerical simulations, the Fuertes model significantly overestimates the wake expansion coefficient $k^{*}$ at $I_u>15\%$. Based on this finding, a logarithmic relationship is proposed and updated as formalized in Equation (37) in this study. The newly developed engineering wake model, termed “logIu”, incorporates this logarithmic relationship, with “$I_u$” denoting the dependence of the model on longitudinal turbulence intensity.

$$k^{*}=\left\{\begin{array}{ll}0.014,\hfill & I_u\le 5\%\hfill \\ {\displaystyle \frac{1}{35}}\ln (I_u)+0.1,\hfill & I_u>5\%\hfill \end{array}\right.$$

(37)

The nonlinear relationship between the wake expansion coefficient and longitudinal turbulence intensity is influenced by the proximity to the ground or sea surface. As the wake expands downstream, it may interact with the ground, producing a compressive effect on the expanding wake. This interaction can increase wind speeds near the ground while decreasing them at hub height, effectively diminishing the wake expansion coefficient. At lower turbulence intensities (<12%), the wake expansion coefficient and range are minimal, resulting in weaker compressive effects from the ground. Thus, the wake expansion coefficient increases linearly with an increase in longitudinal turbulence intensity. However, once turbulence intensity surpasses a critical threshold (12%), the ground’s compressive effect becomes pronounced, inhibiting the linear growth of the wake expansion coefficient with turbulence intensity.

5. Validation and Evaluation of the LogIu Wake Expansion Model

5.1. Wind Tunnel Experiment Validation

Figure 12 compares the measured wind speeds at hub height under neutral conditions in a wind tunnel against predictions from various engineering models. Specifically, Figure 12 delineates the wind speeds along the wake centerline at $TSR=5,6$ for the Dou2019 case. The proposed model shows better agreement with the experiment than other models. The Campagnolo model, which is formulated based on wind tunnel wake measurement data, tends to underestimate wake wind speeds at $x>4D$. Similarly, other wake models also show varying degrees of underestimation, with the Jensen model displaying a notably significant deviation. For identical wake widths, the Jensen model distributes the wind speed reduction uniformly across the wake, leading to an overestimation of speeds near the wake centerline and an underestimation at the edges of the wake width.

The prediction accuracy of the wake speed for WiRE-01 of the proposed model is also higher than those of other models (Figure 13). Relative to the logIu model, the Campagnolo model significantly overestimated the wind speed loss at x = 4D by approximately 37% near the wake centerline. In contrast, other models predict lower wake deficits near the wake center, with the Cheng model showing the most substantial underestimation. In this verification case, the longitudinal turbulence intensity of the hub height is I_u = 7%. However, according to the original method described in the Cheng model’s literature, it is only 5.6%. This discrepancy suggests that the Cheng model underestimates wake deficits when using a lower turbulence intensity. Meanwhile, the predictions by the Fuertes model fall between those of the logIu and Cheng models: the Fuertes model utilizes a larger a_* than the logIu model in wind tunnel scenarios, resulting in a smaller calculated near-wake length and a faster predicted wake recovery, whereas the Cheng model employs a higher turbulence intensity, leading to a quicker wake recovery rate than that of Fuertes.

Figure 14 shows the predicted versus actual wind speeds for various models under different atmospheric conditions. The wake measurement data from the Hancock 2014 study highlights the influence of atmospheric stability on wind turbine wakes: at $x=10D$, the maximum wake deficit is only 15% under unstable conditions, whereas it is 25% under stable conditions. Among the models, the logIu model demonstrates the highest accuracy in predicting wakes. Under non-neutral conditions, the Campagnolo model tends to overestimate wake deficits, failing to deliver accurate predictions at $x=2~5D$. This indicates that despite being optimized with wind tunnel data and performing adequately under neutral conditions, the Campagnolo model lacks satisfactory prediction capability for wind turbine wakes under non-neutral conditions. Similar to neutral conditions, the Jensen, Fuertes, and Cheng models markedly underestimate wake deficits, especially in the near-wake region and under stable atmospheric conditions.

5.2. Field Observation Experiment Validation

5.2.1. Neutral Conditions

Figure 15 evaluates the wake velocity models against measurement data at locations $x=2.5D$, 5.5D, and $8D$ for the Nibe B wind turbine under neutral conditions. The results reveal that the wind speeds predicted by the logIu model align most closely with the actual measurements. The accuracy of the Fuertes and Cheng models is slightly less precise than that of logIu, although the discrepancy is not marked. The Jensen model still underestimates the wake deficits at the wake center ($y\approx 0$) and overestimates it at the wake edge ($y\approx 0.5D$). Because the parameter ${\alpha }_{*}=0.952$ in the Campagnolo model is optimized through wind tunnel tests and is much smaller than the recommended value of 3.6 for actual wind farms, the Campagnolo model predicts a higher value for $x_0$ and thus overestimates wake deficits.

5.2.2. Non-Neutral Conditions

Figure 16 compares the model-predicted wind speeds under different stability conditions with four wind turbine field measurement values. Overall, aside from overestimating the wind speed loss near the wake center under strong unstable conditions ($x=5D$) and strong stable conditions ($x=2D$), the logIu model generally offers reliable predictions across different atmospheric stabilities. The Fuertes and Cheng models exhibit a slightly lower overall predictive accuracy than logIu, with a more pronounced discrepancy evident in the Danwin scenario. The “top hat” shape of the Jensen model still exhibits poor accuracy in predicting wakes in the near-wake region (up to 35% at $2~2.5D$) and overestimates wake deficits at $y=0.5D$. In contrast, the Campagnolo model significantly overestimates wake deficits in these cases, with overestimation ratios reaching up to 170% at $x=3D$ and 30% to 40% at $x=6.1D$, illustrating substantial predictive errors in actual wind farm applications.

5.3. Overall Model Evaluation

Figure 17 compares the RMSE ($U/U_{\infty }$) values for various models across different operational environments. The wake velocities predicted by the logIu and Campagnolo models are closer to the measured values of wind turbine wakes in the wind tunnel compared to the predictions of other engineering models. The RMSE values for logIu and Campagnolo models range from 0.04 to 0.06, which is substantially lower than the 0.07 to 0.09 range observed for alternative models. The prediction accuracy of logIu exhibits a slight decline in actual wind farms, with the corresponding RMSE increasing from approximately 0.04 in the wind tunnel to 0.063. This value remains marginally lower than the 0.07 RMSE associated with Fuertes’s model. The decrease in prediction accuracy may be attributed to the increased complexity of operational environments in real wind farms and potential limitations in measurement instrumentation precision. The prediction accuracies of Fuertes and logIu in actual wind farms are similar because the functional relationships between the wake expansion coefficient and turbulence intensity for the two models demonstrate convergence at $5\%\le I_u\le 12\%$. The prediction accuracy of the Cheng model is slightly lower than that of the Fuertes model, but both are higher than that of the Jensen model considering atmospheric stability. For all cases, logIu has the highest accuracy in wake prediction, followed by the Fuertes and Campagnolo models, with the Jensen model exhibiting the lowest accuracy in wake velocity predictions.

The RMSE ($U/U_{\infty }$) values of various models under different atmospheric stability conditions are shown in Figure 18. The logIu model demonstrates the lowest RMSE across all examined stability conditions, indicating superior overall prediction accuracy. Except for Campagnolo, the reliability of the investigated models in predicting wake velocities generally decreases with increasing stability. In both wind tunnel simulations and real atmospheric boundary layer conditions, increasing stability is associated with a decrease in the boundary layer height. This phenomenon allows wind turbines to exceed the surface layer, where the Reynolds stresses are relatively constant. Above the surface layer, the Reynolds stresses diminish with increasing height, resulting in more pronounced disparities in the turbulence intensity between the upper and lower sections of the wind turbine rotor sweep area. The models studied in this study assume that the relative wind speed loss is symmetrically distributed along the wake centerline, which is not conducive to capturing the increasing difference in turbulence intensity between the upper and lower parts of the rotor sweep area as the stability increases. This limitation in the underlying assumptions of the models may contribute to the observed decrease in predictive accuracy under conditions of heightened atmospheric stability.

6. Conclusions

This study examines the influence of ground effects on wake expansion and proposes a novel nonlinear wake expansion model that incorporates both atmospheric stability and ground effects. The proposed model is calibrated using extensive datasets comprising wake measurements from wind tunnel experiments, field observations, and numerical simulations. A comparative analysis is conducted between the wake wind speeds predicted by the proposed model and those of other typical models. The principal findings of this study are summarized as follows.

(1) Ground effects tend to suppress wake expansion at far-wake locations, and this suppression becomes more pronounced under high turbulence. This ground-wake interaction at far-wake locations, particularly under high turbulence intensity, has the potential to induce a downward displacement of the wake centroid and thus suppress wake expansion at the hub height.

(2) The experimental data and simulation results indicate that the previously assumed linear relationship between the wake expansion coefficient and turbulence intensity is invalid at high turbulence levels. Instead, the wake expansion coefficient exhibits a logarithmic relationship with longitudinal turbulence intensity.

(3) The proposed logIu model demonstrates superior overall accuracy in predicting wake wind speeds, with corresponding RMSE values ranging from 0.04 to 0.063. The prediction accuracy of the wake wind speeds across various models generally exhibits an inverse relationship with increasing atmospheric stability. This trend suggests that wake prediction is the most challenging task under stable atmospheric conditions.

These findings contribute to the advancement of wake modeling techniques and provide valuable insights into the complex interactions among atmospheric conditions, ground effects, and wake behavior in wind energy applications.