Large-eddy simulations (LES) are performed on the flow over a wind farm sited behind an abrupt rough-to-smooth surface roughness jump. The change in surface roughness affects both the first-order and second-order turbulent statistics. The usual deficit, i.e., the difference between the velocities upstream of the entire wind farm and downstream of a turbine, attains negative values close to the ground, which makes it difficult for modeling within the usual Gaussian radial-shape framework. A different definition, i.e., the difference in velocity at the same location with and without a turbine on a heterogeneous surface, is always positive and is amenable to Gaussian shape-based modeling. For the setup considered here, wind farms sited downstream of a surface roughness jump produce more power than a wind farm sited on a homogeneously rough surface. This increase is primarily because of the larger power generated by the downstream turbines and only slightly due to the increased power of the first-row turbine. The farm performance is affected by the distance between the abrupt change in surface roughness and the position of the first row of turbines. The wind farm performance is also dependent on the aerodynamic roughness upstream of the surface roughness jump. Two single-turbine analytical models and three wake-merging strategies are evaluated for their ability to predict the velocity deficits. A corrected form of the standard Gaussian model with a recently proposed wake-merging methodology, applicable for a varying background field, is found to be insensitive to the tunable model parameter and is consistently in line with the LES results.

The flow within the atmospheric boundary layer (ABL) directly impacts the energy extracted by a wind farm. Modern wind farms often span vast horizontal distances, over which the land is likely to undergo a change in its type and use. This leads to spatial variability of the momentum fluxes imposed by the ground surface on the ABL and wind farm flow. The momentum flux imposed by the ground is parameterized by the aerodynamic roughness height (z0). Therefore, it is crucial to investigate the influence of surface roughness heterogeneities on ABL flow and the consequences for wind farms located on such heterogeneously rough terrain.6,27

The impact of surface roughness heterogeneities on ABL flow has been investigated through several experiments or field observations7,9 and numerical simulations,1,27,33,37 with the most commonly studied configuration being that of an abrupt transition between a rougher and smoother surface in the normal direction to that of the mean wind. These investigations have demonstrated the formation of an internal boundary layer (IBL) at the roughness transition point, where the change in surface roughness modulates the mean streamwise velocity and turbulent fields. Several analytical models have been formulated for the mean streamwise wind speed1,9,13,14,28,30 and for the turbulence intensity (TI)27 over surfaces with heterogeneous roughnesses.

Several experimental and numerical studies have been reported of wind turbines/farms and are summarized in the review articles by Stevens and Meneveau38 and Porté-Agel et al.32 A number of simplified analytical models have been developed to predict the evolution of the wake of an isolated wind turbine. Among these, the model by Bastankhah and Porté-Agel2 particularly has been widely adopted because it conserves mass as well as momentum and assumes a Gaussian shape in the radial direction starting at moderate distances downstream of the turbine, which approximates the actual shape very closely. Engineering wind farm wake models (“bottom–up” models following the terminology of Stevens et al.40) rely on strategies to merge the wakes of multiple turbines together and come up with a composite wake. The oldest examples of these wake models that account for interactions between multiple turbines are the linear merging26,29 and quadratic merging.23,47 A new wake merging strategy was proposed recently by Lanzilao and Meyers (LM),25 which, unlike previous models, allows for the background velocity field on which the wind turbine wake-induced deficits are merged, to be non-uniform. Almost all of the above-mentioned numerical simulations or analytical modeling studies have considered wind turbines/farms sited on homogeneously rough surfaces.

Studies involving surface roughness heterogeneities and wind turbines have been limited in number. Kethavath et al.24 studied the evolution of an isolated wind turbine wake with the turbine sited downstream of an abrupt transition in z0. The surface roughness heterogeneity was shown to significantly alter the turbulent statistics. Furthermore, a few corrections to the Gaussian model were proposed, which account for a non-uniform incoming flow and incorporate the effect of a few terms in the momentum deficit transport equation that are usually neglected by other single-turbine models. The model by Kethavath et al.24 was shown to provide better predictions of the large-eddy simulation (LES) results of a single wind turbine wake sited behind an abrupt surface roughness jump than the Gaussian model.

Several studies have investigated coastal effects on offshore wind farms as reviewed in Schulz-Stellenfleth et al.36 For example, the Anholt wind farm in Denmark has been studied using a combination of mesoscale simulations, supervisory control and data acquisition (SCADA) data, and Reynolds-averaged Navier–Stokes (RANS) simulations for understanding the impact of horizontal wind-speed gradients, the effect of the distance between the coast and turbines, and assessing the accuracy of wake models in predicting power losses and annual energy production.31,45 The work of Van Der Laan et al.45 showed that the distance of the leading turbine in a wind farm array to the shoreline affects the wind speed seen and the power generated by the turbine. Recently, Goit and Önder18 investigated the effect of the land–sea transition on the performance of an offshore wind farm. The first row of turbines was placed 25 turbine diameters downstream of the land–sea transition. Incorporating the effect of this surface roughness transition into the simulations led to an enhanced wake recovery, greater inflow turbulence levels, and increased the kinetic energy entrainment. For the specific wind farm layout and surface roughnesses considered, the wind farm placed near the land–sea transition was shown to produce 17% higher power than a wind farm placed further offshore, far from the surface roughness transition.

A key input to wake models for wind farms is the wake width behind each turbine, usually denoted by σ. While analytical modeling of the wake width is a topic of active research,11,12,44 the simplest and most widely used approximation is a linear fit, i.e., σ = σ 0 + k * ( x x t ), where x and xt are the streamwise coordinate and the location of the turbine, respectively, and the fitting parameters are k * (wake growth rate) and σ0 (the initial wake width). Several modeling studies4,21,29,42 have focused on understanding the impact of and developing empirical relations for the wake growth rate, k *, for example, as a function of the turbulence intensity upstream of a turbine. The effect of σ0 has been explored to a limited extent, with values proposed by Bastankhah and Porté-Agel3 being used for most studies. In this paper, we explore the impact of σ0 systematically for wind farms sited on a homogeneously rough surface as well as for wind farms sited downstream of an abrupt change in surface roughness.

In this paper, we study the flow in wind farms sited downstream of a rough-to-smooth abrupt surface roughness transition using large-eddy simulations (LES). Numerical simulations are conducted to study the effect of three parameters, namely, the distance between the surface roughness transition and the first turbine row, the surface roughness values (z01, z02), and their ratio, m = z 01 / z 02. Analytical models for the mean streamwise velocity deficit inside the wind farm are investigated. Specifically, it is of interest to investigate the corrected-Gaussian24 single-turbine model and the wake merging methodology of Lanzilao and Meyers,25 as they are designed to be applicable to the spatially varying background velocity fields that arise on account of surface roughness heterogeneities. We build upon our previous work, which focused on the effect of an abrupt rough-to-smooth surface roughness jump on the ABL14,27 and on the wake of an isolated turbine.24 The present study has similar scope as that of Kethavath et al.24 but for a thirty-turbine wind farm. Additionally, this paper builds upon the work of Goit and Önder,18 which did not consider a systematic variation of parameter values and the evaluation of engineering wake models. Furthermore, we also focus on the effect of σ0 as a key tunable parameter in wind farm wake models and make recommendations for its appropriate values.

The LES methodology, cases simulated, and details of the analytical models evaluated in this study are outlined in Sec. II. Results of the LES are presented and discussed in Sec. III. Evaluation of analytical models to predict the mean velocity deficit in wind farms placed downstream of an abrupt change in surface roughness is discussed in Sec. IV. Summary and conclusions of the work are provided in Sec. V.

Large-eddy simulations (LES) are employed to simulate the flow in wind farms situated over heterogeneously rough surfaces. The incompressible Navier–Stokes equations are solved under the assumption of truly neutral conditions for the atmospheric boundary layer. The influence of the Coriolis force in the atmosphere is not considered in these simulations. A standard concurrent-precursor simulation methodology is utilized, wherein two computational domains are employed.15 The precursor domain has a uniform surface roughness z01 without any turbine sited in it, and the flow here is driven by a constant pressure gradient ( u * 2 / L z). Here, u * is the friction velocity, and Lz is the domain size in the vertical direction. The main domain includes the features of interest, i.e., surface heterogeneity and/or a wind farm. Each domain has the same size ( L x × L y × L z ) along the streamwise (x), spanwise (y), and vertical (z) directions, respectively. The aerodynamic roughness in the main domain is z01 up to x s = 10 d (where d is the diameter of each turbine) and z02 beyond it in cases involving a surface roughness heterogeneity.

The last quarter of the main domain is referred to as the “fringe” region, where an additional forcing term nudges the flow toward the same conditions as in the precursor domain. The aerodynamic surface roughness in the fringe region is the same as in the upstream and precursor domains, i.e., z01. An array of turbines (shown in Fig. 1) are placed in the main domain with three columns and ten rows. The first row of turbines is at a streamwise distance of xt from the inlet plane, and the columns are located at 2d, 6d, and 10d in the spanwise direction. The spacings between the wind turbines along the streamwise and spanwise are sx and sy, respectively. Turbine row n is at x t , n = x t + ( n 1 ) s x, so that x t , 1 is identical to xt. We use spacings s x = 5 d and s y = 4 d, a fixed value of the local thrust coefficient C T = 0.9 (nominal thrust coefficient of C T = 0.6), and hub-height z t = 0.6 d. All these values are reasonable approximations of actual wind farms and have been used in previous idealized studies of wind farm wakes.15,24,39 While not pursued in this paper, sensitivity of the results to these parameters should be studied in future work.

FIG. 1.

Schematic of the main domain (not to scale) presented in the yz and xz planes. Meaning of “HomogT,” “HNT,” and “HTx” is shown in Table I. sx and sy denote the spacing along the streamwise and spanwise directions, respectively. xt and xs are the locations of the first turbine row and surface roughness jump, respectively. z01 and z02 are the two aerodynamic surface roughness values. The last part of the domain is referred to as the “fringe” region, where the surface roughness was set to be z01.

FIG. 1.

Schematic of the main domain (not to scale) presented in the yz and xz planes. Meaning of “HomogT,” “HNT,” and “HTx” is shown in Table I. sx and sy denote the spacing along the streamwise and spanwise directions, respectively. xt and xs are the locations of the first turbine row and surface roughness jump, respectively. z01 and z02 are the two aerodynamic surface roughness values. The last part of the domain is referred to as the “fringe” region, where the surface roughness was set to be z01.

Close modal

Both the horizontal (x, y) directions are assumed to be periodic. The top boundary (z = Lz) has free-slip conditions, while the bottom boundary (z = 0) has a shear stress specified via a wall model. No-penetration vertical velocity conditions are applied to both walls.

The simulations are carried out using the in-house “PadeOps-igrid” code,41 which employs Fourier-spectral discretization in the horizontal (x, y) directions and a sixth-order staggered compact finite-difference scheme in the vertical direction. A third-order Runge–Kutta method is employed for the time advancement. All simulations maintain a constant Courant-Friedrichs-Lewy (CFL) value of 0.8, which is within the stability limit of the time stepping scheme. The anisotropic minimum dissipation (AMD) model35 is used for subgrid-scale closure. A standard actuator-disk model without rotation (ADM-NR)48 and with uniform distribution of the axial thrust force is used to model the turbine forces.22 The viscous sub-layer in this flow is too thin to be resolved by moderately sized grids due to the high Reynolds number ( R e 10 8) based on the boundary layer height, the velocity at the hub-height of a typical wind turbine, and the viscosity of air. Consequently, the region where the direct influence of viscosity is significant is minuscule, and resolving this would require very large grids. This problem is sidestepped by dropping the viscous terms from the Navier–Stokes equations and introducing an additional shear stress at the bottom wall using a wall model. The wall model proposed by Bou-Zeid et al.5 is employed here. Unlike several other wall models, this wall model does not assume homogeneity in the aerodynamic roughness of the underlying surface, is fully local, and hence, is applicable over heterogeneously rough surfaces.5,24,27

The PadeOps-igrid code employed here has been extensively validated and utilized for various studies previously, including those involving LES of stratified and unstratified turbulent boundary layers, wind turbines, and wind farms.15–17,19,20 In our recent previous work,27 we have validated the code for a high-Re half channel flow (as a surrogate for an ABL) over an abrupt rough-to-smooth surface roughness transition by comparing with the wind tunnel experimental data of Chamorro and Porté-Agel8 and demonstrating good agreement for first-order as well as second-order turbulent statistics. Furthermore, in this paper (in  Appendix B), we validate the capability of the code to accurately simulate wind farm wakes by reproducing the experimental data of Chamorro and Porté-Agel10 for a miniature-turbine wind farm sited on homogeneously rough ground surface.

Five sets of LES are conducted, which are presented in Table I. In the first set, illustrated in Fig. 1, case HT1 includes a wind farm with the first turbine at a short distance x t x s = 4 d from the abrupt transition in surface roughness. This case is compared to case HomogT, where the wind farm is retained, but the surface heterogeneity is eliminated. In case HNT (heterogeneous; no turbine), the wind farm is eliminated, but the surface heterogeneity is preserved. The surface roughness values z01 and z02 are taken to match the experiments of Chamorro and Porté-Agel,8 and the case HNT reproduces this experimental dataset (described in detail in Mondal et al.27). The resulting ratio of surface roughnesses, z 01 / z 02 = 83.3, is of the order of magnitude of a transition from smooth land (e.g., grassland) to sea8,43 or from rough suburban landscape to bare soil or snowy surfaces. The computational domain size for all the cases mentioned is ( L x × L y × L z ) =  ( 96 d × 12 d × 4.5 d ). The simulations are carried out on a uniform grid with ( n x × n y × n z ) =  ( 768 × 96 × 108 ) points in each of the precursor and main domains.

TABLE I.

Summary of cases simulated. Location of abrupt change in surface roughness is xs and location of the first turbine row is xt. Aerodynamic roughness height of the upstream, rough, surface is z01 and that of the downstream, smooth, surface is z02. The actual number of grid points per simulation is twice that mentioned here.

Case Description/purpose ( x s , x t ) / d ( L x × L y × L z ) / d ( n x × n y × n z ) z 01 / d m = z 01 / z 02
HomogT  Homogeneous; turbine  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  NA 
HNT  Heterogeneous; no turbine  ( 10 , )  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT1  Heterogeneous; turbine  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT2  Varying xt  (10, 17)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT3  (10, 20)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT4  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT1-m1  Varying z02 and z01  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  25 
HT1-m3  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  125 
HT4-m1  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  25 
HT4-m3  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  125 
HT1- z 01   (10, 14)  96 × 12 × 4.5  768 × 96 × 108  1.125 × 10 3  25 
HT1-g1  Sensitivity to grid and domain size  (10, 14)  96 × 12 × 4.5  256 × 32 × 30  3.75 × 10 3  83.3 
HT1-g2  (10, 14)  96 × 12 × 4.5  512 × 64 × 60  3.75 × 10 3  83.3 
HT1-g3  (10, 14)  96 × 12 × 4.5  640 × 80 × 90  3.75 × 10 3  83.3 
HT1-g5  (10, 14)  96 × 12 × 4.5  1024 × 128 × 120  3.75 × 10 3  83.3 
HT1-g4    (10, 14)  96 × 12 × 6  768 × 96 × 144  3.75 × 10 3  83.3 
HomogE  Validation with experiments10  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  2.0 × 10 4  NA 
HomogT-z0  Upstream TI value matched with HT4  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  8.75 × 10 4  NA 
Case Description/purpose ( x s , x t ) / d ( L x × L y × L z ) / d ( n x × n y × n z ) z 01 / d m = z 01 / z 02
HomogT  Homogeneous; turbine  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  NA 
HNT  Heterogeneous; no turbine  ( 10 , )  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT1  Heterogeneous; turbine  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT2  Varying xt  (10, 17)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT3  (10, 20)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT4  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  83.3 
HT1-m1  Varying z02 and z01  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  25 
HT1-m3  (10, 14)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  125 
HT4-m1  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  25 
HT4-m3  (10, 30)  96 × 12 × 4.5  768 × 96 × 108  3.75 × 10 3  125 
HT1- z 01   (10, 14)  96 × 12 × 4.5  768 × 96 × 108  1.125 × 10 3  25 
HT1-g1  Sensitivity to grid and domain size  (10, 14)  96 × 12 × 4.5  256 × 32 × 30  3.75 × 10 3  83.3 
HT1-g2  (10, 14)  96 × 12 × 4.5  512 × 64 × 60  3.75 × 10 3  83.3 
HT1-g3  (10, 14)  96 × 12 × 4.5  640 × 80 × 90  3.75 × 10 3  83.3 
HT1-g5  (10, 14)  96 × 12 × 4.5  1024 × 128 × 120  3.75 × 10 3  83.3 
HT1-g4    (10, 14)  96 × 12 × 6  768 × 96 × 144  3.75 × 10 3  83.3 
HomogE  Validation with experiments10  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  2.0 × 10 4  NA 
HomogT-z0  Upstream TI value matched with HT4  ( , 5 )  96 × 12 × 4.5  768 × 96 × 108  8.75 × 10 4  NA 

In the second set of simulations, the position of the first row in the wind farm ( x t = 17 d, 20d, and 30d) is changed. These cases enable an assessment of how the position of the turbines relative to the location of the abrupt jump in surface roughness impacts the results.

In the third set, simulations are carried out for different ratios of the surface roughness values, i.e., m = z 01 / z 02. Without changing the upstream surface roughness (z01), simulations are conducted for two additional values of m = 25 and m = 125 (cases HT1-m1 and HT1-m3). In case HT1- z 01 , we keep the z02 same as in case HT1, but the ratio is m = 25.

In the fourth set, the effect of grid size on the results is evaluated by comparing cases HT1-g1, HT1-g2, HT1-g3, and HT1-g5 to the baseline case HT1. This series of simulations demonstrates that 768 × 96 × 108 points are adequate for obtaining converged results as shown in  Appendix A. An additional simulation, case HT1-g4  , is carried out to demonstrate that the vertical domain size is adequately large.

Case HomogE is used to validate our numerical framework for its ability to accurately simulate wake interactions in a wind farm sited on a homogeneously rough ground surface by demonstrating ( Appendix B) good agreement with the experimental results of Chamorro and Porté-Agel.10 

Finally, a comparison is made between cases HomogT, HT4, and HomogT-z0, which has a homogeneously rough surface but with z0 adjusted such that the turbulence intensity upstream of the first turbine matches the value seen in the HT4 case. This comparison is discussed in  Appendix C and makes it clear that the effects of an abrupt surface roughness transition upstream of a wind farm are not explained only by a lower upstream turbulence intensity but are qualitatively different.

All quantities are normalized using the velocity scale u * and the length scale d. Typical values for these may be considered to be u * = 0.45 m/s and d = 100 m. Simulations are conducted for 180 time units, where one time unit equals d / u *. The flow field is given 60 time units to reach a statistically stationary state, and the remaining 120 time units are utilized to calculate statistical averages. These spin up and averaging durations are consistent with prior studies.15,39 The results of the precursor are averaged in both time and the horizontal planes, whereas the statistics for the main domain are averaged in time and along the vertical planes located at turbine mid-span locations ( y = 2 d, 6d, and 10d) of the three turbine columns.

The analytical model proposed by Bastankhah and Porté-Agel2 (called Gaussian model henceforth), which has been widely adopted and successfully used, assumes that the wake velocity deficit decays in the streamwise (x) direction and follows a self-similar Gaussian profile in the radial (y, z) directions. In this paper, we consider two different definitions of velocity deficit. The first is the difference between the incoming upstream velocity and the downstream velocity at any point, i.e., D U ( x , y , z ) = U ( x u p , y , z ) U ( x , y , z ). The second is the difference in the velocity at the same downstream location between case HNT and a heterogeneous case (e.g., HTx), i.e., Δ U ( x , y , z ) = U 0 ( x , z ) U ( x , y , z ). The second definition is necessitated because, as shown below, the first definition attains negative values close to the ground for the heterogeneous cases. The second definition of velocity deficit is always positive for all cases and is considered for analytical modeling purpose in this paper.

The Gaussian model, derived assuming that the incoming wind is uniform in x and z, predicts the velocity deficit due to turbine n at a downstream point (x, y, z) as
Δ U n U 0 , h , u p = [ 1 1 C T 8 ( σ n / d ) 2 ] exp ( r 2 2 σ n 2 ) .
(1)
For flows over heterogeneous surfaces, the incoming wind U0 is a function of x and z, which is not incorporated in the Gaussian model. The varying background velocity field is incorporated into the Gaussian model with slight modifications, similar to what was described in our previous work (Kethavath et al.24). The velocity deficit due to turbine n is given by
Δ U n U 0 , h , u p = Δ U n ¯ U 0 , h , u p + ( U 0 , h U 0 , h , u p 1 ) ,
(2)
where
Δ U n ¯ U 0 , h , u p = [ 1 1 1 F ( x ) ( σ n / d ) 2 [ C T 8 + γ π ] ] exp ( r 2 2 σ n 2 ) .
(3)
In the above-mentioned equations, U 0 , h , u p , r = ( y y t ) 2 + ( z z t ) 2 and σ n ( x ) are the hub-height velocity upstream of the first turbine, the radial distance from the turbine center and the wake width behind turbine n, respectively. U 0 , h = U 0 ( x , z t ) is the velocity at the hub-height in the HNT case, which is obtained from the model described in Ghaisas.14 The function F ( x ) 0.025405 ( x x t , n ) / d 0.060504 is a correction to the Gaussian model proposed by Kethavath et al.24 Another correction proposed by Kethavath et al.,24 to account for additional terms that are missing from the Gaussian model, gives the non-dimensional constant γ 0.01629 87.33 ζ + 7.738 × 10 4 ζ 2, where ζ is obtained from the relation
ζ = z t δ i ( x t , n ) d 2 ln ( δ i ( x t , n ) z t ) z t d / 2 z t + d / 2 1 U 0 2 ( U 0 U 0 x + W 0 U 0 z ) ( x t , n , y t , z ) d z .
Here, W0 is the vertical velocity in the HNT case and is obtained by imposing continuity with the U0 field.24, δi is the internal boundary layer (IBL) height, calculated using the empirical relation13  δ i ( x ) = z 02 [ 0.75 0.03 ln ( z 02 / z 01 ) ] ( ( x x s ) / z 02 ) 0.8. The last term in Eq. (2) is added only up to the top tip-height, i.e., z z t + d / 2.

For a homogeneous surface, the last term is always zero. Furthermore, if we consider the corrections γ and F(x) as zero, Eq. (2) reduces to the Gaussian model [i.e., Eq. (1)]. We refer to Eq. (2) as the Gaussian-corrected (GC) model.

The above models describe the wake of an isolated wind turbine. Multiple wakes must be accounted for to predict the cumulative wake losses in wind farms. In this study, we consider three wake superposition techniques. Two of the most commonly used wake merging techniques are the linear sum of velocity deficits,26,29
Δ U ( x , y , z ) = n = 1 N t Δ U n ( x , y , z ) ,
(4)
and the square-root of summation of the squares of velocity deficits,23,47
Δ U ( x , y , z ) = [ n = 1 N t Δ U n 2 ( x , y , z ) ] 1 / 2 .
(5)
These techniques have also been termed as linear and quadratic wake merging, respectively. A new wake superposition method developed by Lanzilao and Meyers25 is also considered here, specifically because it accounts for a heterogeneous background velocity field and is expected to be applicable to flows involving surface roughness heterogeneities. The composite wake is given by
Δ U ( x , y , z ) = U 0 ( x , z ) [ 1 n = 1 N t ( 1 Δ U n ( x , y , z ) / U n 1 ( x , z ) ) ] .
(6)
We combine the three wake merging methodologies, namely, “linear” (shortened to L), “quadratic” (Q), and “Lanzilao and Meyers” (LM) with the Gaussian model (shortened to G) to obtain the combinations GL, GQ, and GLM, respectively. To combine the new single-turbine models [Eqs. (2) and (3)] proposed here, we apply the merging techniques to Δ U n ¯ rather than to Δ U n, to get the cumulative Δ U ¯. In other words, Δ U ¯ is given by Eqs. (4)–(6) with Δ U n ¯ replacing Δ U n. Finally, this cumulative deficit is corrected only inside the rotor-region, similar to Eq. (2),
Δ U ( x , y , z ) = Δ U ¯ ( x , y , z ) + ( U 0 , h ( x ) U 0 , h , u p ) .
(7)
These three combinations are termed as GCL, GCQ, and GCLM, where “GC” indicates “Gaussian-corrected.”

Cases HomogT, HNT, and HT1 are discussed first. This allows a comparison of the performance of a wind farm sited behind a surface roughness jump with that of a wind farm sited on homogeneously rough ground surface. The results of the precursor simulations and the profiles upstream of the first-row turbines from the main domain (i.e., x / d = x t , 1 1) for the cases HomogT, HNT, and HT1 are shown in Fig. 2. In the HNT case, where there is no wind farm and hence xt is not defined, the profiles are plotted at the same location as in the HT1 case. Profiles from the left to right are the time-averaged streamwise velocity (U), streamwise turbulence intensity (TI), and negative total shear stress (−TSS). The streamwise turbulence intensity is defined as T I = U U ¯ / U h, where Uh is the velocity at the hub-height 1d upstream of the first turbine row. The total shear stress is the sum of resolved, subgrid-scale and wall model components, TSS = U W ¯ + τ 13 sgs ¯ + τ 13 w m ¯. The profiles of U, TI, and TSS vary by 0.3%, 1.7%, and 5.9%, respectively, between the precursor and the main domain for the HomogT case. This close agreement is expected since the aerodynamic surface roughness is the same ( z 01 ) in the precursor domain and in the main domain in the HomogT case. The small differences can be attributed to induction effects caused by the wind farm. For the cases HT1 and HNT, close to the lower wall, there is an increase in U and a decrease in both TI and TSS compared to the HomogT case. This is because the surface roughness has changed from z01 to z02 at the location 1d upstream of the first turbine row. The profiles of U, TI, and TSS match closely with the precursor and HomogT profiles above a certain height, as is also expected, highlighting the well-controlled nature of the inflow conditions applied in our numerical simulations.

FIG. 2.

Vertical profiles of normalized time-averaged streamwise velocity (left panel), streamwise turbulence intensity (middle panel), and negative of total shear stress (right panel) at 1d upstream of the first-row turbine for cases HomogT, HNT, and HT1 as shown in Table I. “Prec” refers to results of the precursor simulation.

FIG. 2.

Vertical profiles of normalized time-averaged streamwise velocity (left panel), streamwise turbulence intensity (middle panel), and negative of total shear stress (right panel) at 1d upstream of the first-row turbine for cases HomogT, HNT, and HT1 as shown in Table I. “Prec” refers to results of the precursor simulation.

Close modal

1. Mean velocity

Contour of the instantaneous streamwise velocity at both x–z and x–y planes for the case HT1 are shown in Fig. 3. The surface roughness change is marked by the dashed vertical line, and the location of the wind turbines is marked by the solid white line. As a result, wind turbine wakes are seen and identified by the low-speed regions downstream of each turbine. Additionally, vertical profiles of the time-averaged streamwise velocity (U) at mid-span and averaged over the three columns of turbines obtained from the LES are shown in Fig. 4. The results of cases HomogT, HNT, and HT1 are shown at different locations downstream of three rows (i.e., row 1, row 5, and row 9). In the cases with wind turbines, HomogT and HT1, a low-speed region downstream of each turbine is seen, which indicates the turbine wake. In cases HNT and HT1, the surface roughness value changes from z01 (rough) to z02 (smooth) at x s = 10 d. After this step change, an increased streamwise velocity near the bottom wall is seen at all downstream locations compared to case HomogT. As seen in Fig. 4, close to the bottom wall, a slight reduction in the streamwise velocity is seen for case HT1 as compared to case HNT.

FIG. 3.

Contours of the instantaneous streamwise velocity for the HT1 case in the (a) xz plane at the mid-span location and (b) xy plane at the turbine hub-height. The turbine rotor disks are marked by solid white lines. The location of the step change in surface roughness is denoted by the dashed black line.

FIG. 3.

Contours of the instantaneous streamwise velocity for the HT1 case in the (a) xz plane at the mid-span location and (b) xy plane at the turbine hub-height. The turbine rotor disks are marked by solid white lines. The location of the step change in surface roughness is denoted by the dashed black line.

Close modal
FIG. 4.

Vertical profiles of normalized time-averaged streamwise velocity for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The shaded region represents the rotor region.

FIG. 4.

Vertical profiles of normalized time-averaged streamwise velocity for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The shaded region represents the rotor region.

Close modal

2. Turbulence intensity

Vertical profiles of the streamwise turbulence intensity (TI) at the turbine center plane, averaged over the three turbine columns, are shown in Fig. 5. It is defined as the ratio of the standard deviation of the streamwise velocity fluctuations to the mean velocity at the turbine hub-height (i.e., T I = U U ¯ / U h). The results are shown at four downwind locations ( x x t , n ) / d = 1–4 behind the same three rows as previously, i.e., row 1, row 5, and row 9. Under both, HomogT and HT1 conditions, a significant increase in turbulence intensity is observed behind each turbine and between z = z t ± d / 2. At each downstream location, the maximum of the TI profile occurs at the top tip height because of the strong shear and correspondingly large turbulent kinetic energy (TKE) production at that location.

FIG. 5.

Vertical profiles of streamwise turbulence intensity (TI) for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The gray shaded region represents the rotor region.

FIG. 5.

Vertical profiles of streamwise turbulence intensity (TI) for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The gray shaded region represents the rotor region.

Close modal

In cases HNT and HT1, the impact of the surface roughness transition is evident, especially close to the lower wall. The TI values are lower in these cases compared to HomogT, where a higher surface roughness exists throughout the domain. The TI profiles obtained from cases HomogT and HT1 are identical behind row 1 except very close to the bottom wall. Behind rows 5 and 9, the turbulence profiles are appreciably smaller for case HT1 than for case HomogT over a larger portion of the domain in the vertical direction. Thus, it is clear that the effect of reduced surface roughness is limited to the region close to the ground surface ( z / d 0.8) at smaller downstream distances from the step change, but it increasingly affects the entire boundary layer at larger downstream distances.

3. Total shear stress

Vertical profiles of the normalized total negative shear stress (TSS), averaged over the three turbine column midplanes, obtained from the LES, are shown in Fig. 6. The total shear stress is given by the sum of the resolved, subgrid-scale and wall model contributions, TSS = ( U W ¯ + τ 13 sgs ¯ + τ 13 w m ¯ ). The surface roughness transition from z01 to a lower z02 in cases HNT and HT1 results in reduced values of TSS near the bottom wall. In cases HomogT and HT1, the effects of turbines are clearly seen, with an increased −TSS above the turbine hub-height and reduced −TSS below the hub-height in the wake region (roughly z t ± d / 2). Consistent with the profiles of TI, the shear stress magnitudes are smaller for the HT1 case than for the HomogT case. These effects are confined to only close to the bottom wall ( z / d 0.8) behind row 1, to roughly z / d 2.5 behind row 5 and to about z / d 3.5 behind row 9.

FIG. 6.

Vertical profiles of negative total shear stress (TSS) for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The gray shaded region represents the rotor region.

FIG. 6.

Vertical profiles of negative total shear stress (TSS) for cases HomogT, HNT, and HT1 as shown in Table I. Results are shown for four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9). The gray shaded region represents the rotor region.

Close modal

4. Internal boundary layer height

The profiles of TI and TSS shown earlier deviate from each other among the cases HomogT, HNT, and HT1 and from the precursor plots at a certain height. These observations indicate the formation of an internal boundary layer that progressively increases in height with increasing streamwise distance. The evolution of the IBL height is shown for the HomogT, HNT, and HT1 cases in Fig. 7. The IBL height at each x location is calculated as the highest position where the upstream TSS profile and the TSS profile at x differ by more than 10%. Calculating the IBL height based on other profiles, e.g., TI, gives similar values.27 The x coordinate for the HomogT case is translated such that the position of the first turbine row matches between the HomogT and HT1 cases.

FIG. 7.

Formation of internal boundary layer for cases HomogT, HNT, and HT1 as shown in Table I. Here, x adj = x x s for the cases HNT, HT1, and x adj = x 1 d for the case HomogT.

FIG. 7.

Formation of internal boundary layer for cases HomogT, HNT, and HT1 as shown in Table I. Here, x adj = x x s for the cases HNT, HT1, and x adj = x 1 d for the case HomogT.

Close modal

The height of the IBL formed due to the surface roughness jump in the HNT case is zero upstream of the step change and increases steadily with increasing x adj = x x s. The height of the IBL formed due to the wind farm in the HomogT case is zero until the turbine location, jumps abruptly to δ i = 1.2 d at the turbine location, and then increases gradually. For the HT1 case, the IBL is formed due to the surface heterogeneity at xadj = 0 and due to the wind farm starting from x adj = 4 d. The IBL height follows the curve for the HNT case until the first turbine row location and closely follows the HomogT curve after the first row of turbines is encountered. These results indicate that the IBL height due to the wind farm is not affected by the presence of the roughness heterogeneity of the underlying surface.

5. Velocity deficit

In Fig. 8, a comparison of three velocity deficit profiles, evaluated from the time-averaged LES results, is presented. The conventional definition of velocity deficit at any point downstream of a turbine is the difference between the incoming upstream velocity and the velocity at that point, i.e., D U ( x , y , z ) = U ( x u p , y , z ) U ( x , y , z ). Here, we take xup to be equal to x t 1 d. In case HomogT, close to the bottom wall, DU has positive values at all streamwise locations. The velocity deficit peaks at the hub-height and decreases gradually to zero above and below this height. In contrast, for case HT1, the DU profile shows significant negative values close to the lower boundary. This is due to the acceleration of the flow downstream of x s = 10 d due to the surface roughness change, which results in U ( x , y , z ) being larger than U ( x u p , y , z ).

FIG. 8.

Vertical profiles of time-averaged velocity deficit at four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9) for cases HomogT, HNT, and HT1 as shown in Table I. D U ( x , y , z ) = U ( x u p , y , z ) U ( x , y , z ) is the velocity deficit as per the usual definition. Δ U ( x , y , z ) = U 0 ( x , z ) U ( x , y , z ) is an alternate definition of velocity deficit.

FIG. 8.

Vertical profiles of time-averaged velocity deficit at four downstream locations of each turbine for three rows (i.e., rows 1, 5, and 9) for cases HomogT, HNT, and HT1 as shown in Table I. D U ( x , y , z ) = U ( x u p , y , z ) U ( x , y , z ) is the velocity deficit as per the usual definition. Δ U ( x , y , z ) = U 0 ( x , z ) U ( x , y , z ) is an alternate definition of velocity deficit.

Close modal

An alternate definition of velocity deficit is considered, defined as the difference between the velocity at the same downstream locations in cases HNT and HT1, i.e., Δ U ( x , y , z ) = U 0 ( x , z ) U ( x , y , z ). Here, U0 refers to the velocity field obtained from the simulation with heterogeneity but no turbine (case HNT). This velocity field is spatially averaged in the spanwise (y) direction due to homogeneity in that direction. The profiles of Δ U shown in Fig. 8 exhibit non-negative values near the bottom wall and are similar to the DU profiles obtained from the HomogT case. The advantage of the alternate definition of velocity deficit is that, in the region 0 z 2 d that is directly affected by the turbine wakes, it has positive values, unlike the DU profile. As a result, unlike the DU profile, the Δ U profile can be modeled by a function that is always positive (e.g., a Gaussian).

The effect of the distance between the surface roughness jump (xs) and the first turbine row in a wind farm (xt) is studied in this section. The qualitative flow features outlined for x t = 14 d (Case HT1) in Subsection III A are seen in cases HT2, HT3, and HT4 as well, where the distance xt is progressively increased to 17d, 20d, and 30d, respectively. These results are not repeated here for brevity.

The effect of increasing xt on the IBL height is shown in Fig. 9. The IBL height increases gradually from the step change location ( x adj = x x s = 0) until the position of the first turbine row, as seen in Fig. 9(a). Once the first turbine row is encountered, the IBL height abruptly jumps to δ i = 1.2 d for all cases. The magnitude of this abrupt jump reduces across the HT1 through HT4 cases because the IBL upstream of xt, due to the surface heterogeneity alone, grows to a larger height as xt increases. Figure 9(b) shows that the IBL heights downstream of the first turbine row are almost the same for all cases, which indicates that they are not affected by the underlying surface roughness.

FIG. 9.

Internal boundary layer height for the various cases indicated in Table I plotted against (a) x adj / d and (b) ( x x t ) / d. For panel (a), the horizontal axis is x adj = x x s for all HTx cases and x adj = x 1 d for case HomogT. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

FIG. 9.

Internal boundary layer height for the various cases indicated in Table I plotted against (a) x adj / d and (b) ( x x t ) / d. For panel (a), the horizontal axis is x adj = x x s for all HTx cases and x adj = x 1 d for case HomogT. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Close modal

The effect of increasing xt is quantified in terms of the power in Fig. 10. Because all turbines have the same rotor diameter and the same local thrust coefficient, the power of each turbine is calculated as P n = U n 3 ¯, where Un represents streamwise velocity averaged over the rotor disk of turbine n. The overhead bar denotes averaging over time and turbine columns, and the subscript n represents the turbine row within the wind farm. While this power is normalized by u * 3, we refer to it as “normalized absolute” power to distinguish it from the relative power. The relative power for turbine row n is computed by normalizing by the power of the first-row turbine, as P n / P 1.

FIG. 10.

(a) Normalized absolute power and (b) relative power, i.e., power normalized by the power of the first-row turbine (P1), as a function of turbine row number for the HomogT and the different HTx cases. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

FIG. 10.

(a) Normalized absolute power and (b) relative power, i.e., power normalized by the power of the first-row turbine (P1), as a function of turbine row number for the HomogT and the different HTx cases. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Close modal

The normalized absolute power [Fig. 10(a)] of first-row turbines increases with increasing xt because of the acceleration of the flow due to the step change in surface roughness leading to an increased disk-averaged velocity. The first-row turbine powers are 0.7%, 3.5%, 5.5%, and 13.4% higher in the cases HT1, HT2, HT3, and HT4, respectively, than in the HomogT case. In the cases HT1 through HT4, the powers generated by the rest of the turbines (rows 2 through 10) are larger compared to case HomogT. The powers P2 through P10 are slightly different between the HT cases with different xt. This indicates a complicated interaction between turbine wakes and the acceleration of the mean flow due to the surface roughness jump.

The relative powers of the downstream turbines (rows 2 through 10), shown in Fig. 10(b), are smallest for case HomogT. Changing to a smooth surface reduces the wake losses and increases the relative powers of the downstream turbines. As xt increases, the relative powers start reducing, which suggests that as the distance from the step change to the first row becomes large, the wake losses again approach the same level as for a homogeneous surface. This is quantitatively shown in Fig. 11, where the average of the relative powers of the downstream turbines is plotted for different cases. Compared to the HomogT case, the relative powers in the HTx cases are larger by about 19.3% (for HT1) to about 10.8% (for HT4). These results show that for a wind farm sited on a heterogeneous surface, the distance between the surface roughness transition and the first turbine row ( x t x s) is a crucial parameter, and that an optimal value of ( x t x s ) exists, for which the wake losses are minimized.

FIG. 11.

Powers of the downstream turbines averaged over rows 2 10 and normalized by the power of the first-row turbine (P1) for the HomogT and the different HTx cases. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

FIG. 11.

Powers of the downstream turbines averaged over rows 2 10 and normalized by the power of the first-row turbine (P1) for the HomogT and the different HTx cases. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Close modal

The normalized absolute and relative powers for different m = z 01 / z 02 values with the first turbine row at xt = 14 (corresponding to case HT1) are shown in Fig. 12. First, we consider a fixed z 01 = 3.75 × 10 3 d and varying z02 to obtain different m ratios. The values of z02 are 1.5 × 10 4 d , 4.5 × 10 5 d, and 3 × 10 5 d for the m ratios 25, 83.3, and 125, respectively. The normalized absolute and relative powers of the HT1-m1 and HT1-m3 cases are almost the same as for the case HT1. Specifically, Fig. 12(a) demonstrates that the normalized absolute power of the front-row turbine is nearly equal under all conditions. The same is found to be true between cases HT4-m1, HT4, and HT4-m3, shown in Fig. 13. Combined with Fig. 10(a), it is clear that, within the parameter ranges studied, the power of the front-row turbine and the wake losses in the wind farm are dependent on the distance from the heterogeneity, but not strongly dependent on the strength of the heterogeneous roughness jump, as long as the upstream surface roughness is kept the same. It is possible that this conclusion does not hold for much longer (more than ten rows) extended wind farms, which can be studied in the future.

FIG. 12.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), as a function of the turbine row number for different m values. All cases have a fixed z 01 = 3.75 × 10 3 d. The values of m for cases HT1-m1, HT1, and HT1-m3 are 25, 83.3, and 125, respectively.

FIG. 12.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), as a function of the turbine row number for different m values. All cases have a fixed z 01 = 3.75 × 10 3 d. The values of m for cases HT1-m1, HT1, and HT1-m3 are 25, 83.3, and 125, respectively.

Close modal
FIG. 13.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), as a function of the turbine row number for different m values. All cases have a fixed z 01 = 3.75 × 10 3 d. The values of m for cases HT4-m1, HT4, and HT4-m3 are 25, 83.3, and 125, respectively.

FIG. 13.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), as a function of the turbine row number for different m values. All cases have a fixed z 01 = 3.75 × 10 3 d. The values of m for cases HT4-m1, HT4, and HT4-m3 are 25, 83.3, and 125, respectively.

Close modal

Figure 14 shows the results of case HT1- z 01 , where the ratio m = 25 is obtained by keeping z02 unchanged from case HT1 but by changing z 01 = 1.125 × 10 3 d. For case HT1- z 01 , the front-row turbine power is much larger than that in cases HT1 and HT1-m1, which is directly related to the reduced upstream roughness z01. Figure 14(b) shows that a significant difference exists in the relative powers of the downstream turbines for case HT1- z 01 compared to HT1 and HT1-m1.

FIG. 14.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), are plotted as a function of the turbine row number for different m values obtained by changing z01 and z02. Cases HT1 and HT1-m1 have a fixed z 01 / d = 3.75 × 10 3, and m = 83.3 and 25, respectively. Cases HT1 and HT1- z 01 have the same z 02 = 4.5 × 10 5, with m = 83.3 and 25, respectively.

FIG. 14.

(a) Normalized absolute power and (b) relative power, normalized by the first-row turbine power (P1), are plotted as a function of the turbine row number for different m values obtained by changing z01 and z02. Cases HT1 and HT1-m1 have a fixed z 01 / d = 3.75 × 10 3, and m = 83.3 and 25, respectively. Cases HT1 and HT1- z 01 have the same z 02 = 4.5 × 10 5, with m = 83.3 and 25, respectively.

Close modal

To summarize, the surface roughness m affects the performance of a wind farm if the upstream roughness (z01) is changed keeping the downstream roughness (z02) fixed but not if the downstream roughness is changed keeping the upstream roughness fixed.

The ability of different analytical models to predict the flow features described earlier is discussed in this section. We focus on describing the HTx cases for brevity. The other cases are predicted with a similar accuracy and are not discussed here.

The performance of the front-row turbines is considered separately from that of the subsequent, downstream, rows. In the HTx cases, the power generated by the most upstream turbine increases as the distance from the surface roughness jump is increased, as shown in Fig. 10(a). This increase is directly linked to the increased wind speed at the hub-height due to the abrupt change in surface roughness. The velocity field over a heterogeneous surface in the absence of turbines, U 0 ( x , z ), is predicted using the model of Ghaisas.14 This predicted velocity field agrees well with the LES results, as seen in Fig. 15.

FIG. 15.

Disk-averaged mean velocities upstream of the first turbine row in the HTx cases normalized by the disk-averaged mean velocity in the HomogT case obtained from the LES and from the model of Ghaisas,14 with parameters α = 0.001 and β = 0.005. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

FIG. 15.

Disk-averaged mean velocities upstream of the first turbine row in the HTx cases normalized by the disk-averaged mean velocity in the HomogT case obtained from the LES and from the model of Ghaisas,14 with parameters α = 0.001 and β = 0.005. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Close modal

Modeling the wake effects of upstream turbines is critical to accurately modeling the performance of the downstream turbines. A crucial input to wake models is the wake growth rate behind each turbine row.

The evolution of the wake widths along the streamwise direction behind each turbine row for the HomogT and HT1 cases is shown in Fig. 16. At each x, the width is obtained as the geometric mean ( σ n = σ z σ y) of the wake widths along the spanwise (σy) and vertical (σz) directions. The σy and σz are computed from the velocity deficit profiles along y and z, respectively, using the expressions provided in Bastankhah and Porté-Agel.3 

FIG. 16.

Evolution of the width of the wake behind each of the turbine rows in the (a) HomogT and (b) HT1 cases. Symbols are the LES results, and solid lines are linear fits to the data in ( x x t , n ) [1d,4d].

FIG. 16.

Evolution of the width of the wake behind each of the turbine rows in the (a) HomogT and (b) HT1 cases. Symbols are the LES results, and solid lines are linear fits to the data in ( x x t , n ) [1d,4d].

Close modal

Figure 16 shows that the wake width evolution is similar for the HomogT and the HT1 case. In each case, the first row turbine wake width is the smallest, and the wake widths increase monotonically from row 1 to row 10. This is expected since the turbulence intensity is the smallest upstream of the first row and increases as one proceeds into the wind farm. Behind each turbine row, the wake has an almost linear expansion, similar to previous observations.2,3,15 A linear fit, σ n σ 0 , n + k n * ( x x t , n ), gives the growth rate, k n *, and the intercept (initial wake width), σ 0 , n, for each n. The wake growth rates extracted from the LES for the HomogT and the HTx cases are shown in Fig. 17(a). It is clear that the behavior is qualitatively similar for the HomogT and all the HTx cases with small quantitative differences between the HTx cases.

FIG. 17.

Wake growth rates behind different turbine rows in the HomogT and the different HTx cases, plotted as a function of (a) the turbine row number and (b) the added streamwise turbulence intensity, Δ T I. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

FIG. 17.

Wake growth rates behind different turbine rows in the HomogT and the different HTx cases, plotted as a function of (a) the turbine row number and (b) the added streamwise turbulence intensity, Δ T I. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Close modal

Figure 16 further shows that, for a given turbine row, the σn values are larger for the HT1 case than for the HomogT case. The reason for this difference can be traced to the fact that the deficit Δ U for the HT1 case is larger than the deficit DU for the HomogT case, particularly close to the ground below the turbine hub-height, as seen in Fig. 8. The difference in the background velocity, namely, the incoming upstream velocity profile for HomogT and the spatially evolving profile from the HNT case for the HT1 case, lead to the deficits DU and Δ U being different, which in turn leads to differences in the σn and k * values between the HNT and HT1 cases. The same observations hold for the other HTx cases and explain the larger k * values seen in Fig. 17(a) for the HTx cases over the HomogT case.

The wake growth rate is usually found to be correlated with the streamwise turbulence intensity.15,29 This does not hold directly for HTx cases since the k * values are larger [Fig. 17(a)], whereas the TI values are smaller (Fig. 5) in the HTx cases compared to the HomogT case, particularly for downstream rows. The wake growth rate is seen in Fig. 17(b) to correlate well (correlation coefficient of 0.9) with the added turbulence intensity, Δ T I = T I 2 T I 0 2. Here, TI0 is the turbulence intensity in the HNT case and, in practical applications, can be obtained as a function of x and z from the model described by Mondal et al.27 for cases involving an abrupt surface roughness jump. For wind farms on homogeneously rough surfaces, TI0 profile is the same at all streamwise locations, equal to the profile upstream of the wind farm, and the Δ T I is the same as the added turbulence intensity used in previous studies.4,29

Two single-turbine wake models (Gaussian and Gaussian-corrected) coupled with three wake merging methods (linear, quadratic, and Lanzilao and Meyers), described in Sec. II C, are compared with the LES data. The normalized velocity deficits at the turbine hub-height and in the turbine center-plane along the spanwise direction, i.e., Δ U ( x , y t , z t ) / U 0 , h , u p for the HomogT and HT1 cases, are shown in Figs. 18 and 19, respectively. To quantitatively assess the accuracy of the different methods, the mean absolute error between the LES and the model predictions is calculated and shown in Figs. 18(d)–18(f) and 19(d)–19(f). The k * values are those extracted directly from the LES, while different choices are considered for σ0. The lines labeled as “Act σ0” use the actual values extracted from the LES. The lines labeled as “Fixed σ0” use a uniform, fixed, value for each turbine, varying between 0.34 and 0.44 in steps of 0.01. The large thickness of these lines captures the upper and lower bounds of the velocity deficit predictions as σ0 is varied over this range.

FIG. 18.

Evaluation of analytical models for the HomogT case. Normalized velocity deficits at the turbine hub-height and in the turbine mid-span plane with (a) linear, (b) quadratic, and (c) Lanzilao and Meyers wake-merging methodology. Mean absolute errors between the LES and model results for (d) linear, (e) quadratic, and (f) Lanzilao and Meyers merging.

FIG. 18.

Evaluation of analytical models for the HomogT case. Normalized velocity deficits at the turbine hub-height and in the turbine mid-span plane with (a) linear, (b) quadratic, and (c) Lanzilao and Meyers wake-merging methodology. Mean absolute errors between the LES and model results for (d) linear, (e) quadratic, and (f) Lanzilao and Meyers merging.

Close modal
FIG. 19.

Evaluation of analytical models for the HT1 case. Normalized velocity deficits at the turbine hub-height and in the turbine mid-span plane with (a) linear, (b) quadratic, and (c) Lanzilao and Meyers wake-merging methodology. Mean absolute errors between the LES and model results for (d) linear, (e) quadratic, and (f) Lanzilao and Meyers merging.

FIG. 19.

Evaluation of analytical models for the HT1 case. Normalized velocity deficits at the turbine hub-height and in the turbine mid-span plane with (a) linear, (b) quadratic, and (c) Lanzilao and Meyers wake-merging methodology. Mean absolute errors between the LES and model results for (d) linear, (e) quadratic, and (f) Lanzilao and Meyers merging.

Close modal

Several observations can be made from Figs. 18 and 19. First, for both (HomogT and HT1) cases, linear wake merging generally overestimates the velocity deficits [Figs. 18(a) and 19(a)]. On the other hand, quadratic wake merging generally underestimates the deficits [Figs. 18(b) and 19(b)]. The Lanzilao and Meyers wake merging methodology gives the closest predictions for both cases, as shown in Figs. 18(c) and 19(c).

The results of linear and quadratic wake merging methodologies are very sensitive to σ0, when the same value is used for each turbine wake. In the HomogT case, the errors of the GL and GCL models vary from roughly 7% to 19% and from 5.5 % to 10%, respectively [Figs. 18(d) and 19(d)]. In the HT1 case, the errors similarly vary from around 5% to about 12%. Considering the GQ and GCQ models, the errors in predictions of the HomogT [Fig. 18(e)] and HT1 [Fig. 19(e)] cases reach as high as 20% and 16%, respectively. The Lanzilao and Meyers merging is not as sensitive to σ0 as seen quantitatively [Figs. 18(f) and 19(f)], with the errors varying between roughly 5% and 10% for the HomogT and HT1 cases.

The sensitivity of the errors to σ0 is primarily because of the sensitivity of the predictions of the first few (up to 5) turbines in each wind farm. By the seventh row ( x / d = 35), the model predictions are almost the same regardless of σ0. This is seen from the vertical spread in the green and blue solid lines in Figs. 18(a)–18(c) and 19(a)–19(c).

In most cases, the Gaussian-corrected wake model gives a smaller error and closer prediction of the LES results than the corresponding Gaussian wake model. This is particularly true in the HomogT case for linear and Lanzilao and Meyers merging [Figs. 18(d) and 18(f)] and in the HT1 case for the quadratic and Lanzilao and Meyers merging [Figs. 18(e) and 18(f)].

Carefully evaluating all model combinations, it is observed that the GCLM method with a fixed, uniform, σ0 set to any value between 0.38 and 0.42 consistently gives small errors of around 5%–6% across all cases. Thus, it is possible to make a clear recommendation of using the corrected Gaussian model with the Lanzilao and Meyers merging method and a fixed uniform value of σ0. The mean absolute errors for the other three HTx cases are compiled in Table II. It is clear that the recommended method gives accurate predictions, with errors being consistently around 5%–6%.

TABLE II.

Relative errors (expressed as percentages) between the LES results and predictions of different model combinations for cases HT2, HT3, and HT4, with a fixed σ0 and k * taken from the LES. “G” and GC refer to the Gaussian and Gaussian-corrected single-turbine models, respectively. The rows refer to the three wake-merging techniques used in the wind farm model. Cases HT1, HT2, HT3, and HT4 have ( x t x s ) / d = 4 , 7 , 10, and 20, respectively, and a fixed surface roughness ratio of m = 83.3.

Fixed σ 0 = 0.42 HT2 HT3 HT4
G GC G GC G GC
Linear  4.8  5.8  4.9  6.5  5.1  9.9 
Quadratic  17.6  10  17.6  10.1  17.4  7.9 
LM  6.4  5.1  6.3  5.4  6.2  6.3 
Fixed σ 0 = 0.42 HT2 HT3 HT4
G GC G GC G GC
Linear  4.8  5.8  4.9  6.5  5.1  9.9 
Quadratic  17.6  10  17.6  10.1  17.4  7.9 
LM  6.4  5.1  6.3  5.4  6.2  6.3 

The evolution of the wake of a wind farm sited on a surface undergoing an abrupt change in aerodynamic roughness from a relatively high value to a relatively low value is studied using large-eddy simulations. The flow in the heterogeneous case is compared to the flow over a homogeneous surface with a wind farm (HomogT) on the one hand and to the flow over a heterogeneous surface without wind farm (HNT) on the other. The change in surface roughness leads to significant differences in the different turbulent quantities. Specifically, the mean velocity is larger, while the magnitudes of the turbulence intensity and the vertical shear stress are smaller in wind farms sited downstream of an abrupt rough-to-smooth roughness transition as compared to a wind farm sited on a homogeneously rough surface. The vertical extent over which these differences are significant increases with increasing downstream distance from the roughness jump. An internal boundary layer develops starting from the position of surface roughness jump and affects the turbulent statistics in the wake of the wind farm.

The velocity deficit (denoted DU), when defined in the usual manner as the difference between the velocity upstream of and inside the wind farm, attains negative values close to the ground for a wind farm sited downstream of an abrupt change in surface roughness. An alternate definition (denoted Δ U), namely, the difference between velocities on a heterogeneously rough surface in the absence of and in the presence of the wind farm, gives positive values in the rotor region. This alternate definition of velocity deficit is more beneficial since it can be modeled as, e.g., a Gaussian in the radial direction, which is always positive. The Δ U values are found to be larger than the DU values due to differences in the “background” velocity field. As a result, the wake width values (σ) and the wake growth rates ( k *) calculated from wind farms sited on homogeneously rough surfaces and surfaces with a roughness jump are different, with the values being larger in the latter cases.

The position of the wind farm with respect to the step change in surface roughness has a significant impact on the wake losses and the power generated. The power generated by the front-row turbines changes by 0.7%–13.4%, while the normalized powers of the downstream turbines, which are indicative of the wake losses, change by about 10.8%–19.3% over the cases simulated. Placing the first row of the wind farm 4 diameters downstream of the surface roughness jump is found to be optimal. As the distance to the first turbine row increases, the wake losses of the downstream turbines systematically increase. For very large distances, the wake losses of the downstream turbines approach the same level as in the HomogT case.

Over the range of parameters simulated, changing the aerodynamic roughness height of the downstream (smooth) surface does not affect the wake losses. On the other hand, changing the roughness of the upstream (rough) surface affects the relative power production by the turbines and the wake losses.

Finally, a detailed description of analytical modeling considerations for wind farms sited on homogeneously rough surfaces and on surfaces with an abrupt roughness transition is presented. First, the power of the upstream turbines for different cases is directly related to the velocity seen by the first row turbines. This velocity is shown to be predicted well by our previously developed model (Ghaisas14). Second, the wind turbine wake growth rates in wind farms sited behind a surface roughness jump are shown to be correlated with the added turbulence intensity, T I 2 T I 0 2. Here, TI0 refers to the turbulence intensity of the background flow due to only the step change in surface roughness (the HNT case) and can be prescribed by the simple model described in Mondal et al.27 

Third, different analytical modeling choices are evaluated for their ability to predict the wake losses in wind farms sited downstream of an abrupt surface roughness jump. Two models for the wake of an isolated turbine are coupled with three distinct wake superposition methods. Results of linear and quadratic wake merging strategies are very sensitive to σ0, which can be considered to be a tunable parameter in the model. The Lanzilao and Meyers wake merging strategy is less sensitive to σ0. A value between 0.38 and 0.42 is found to be suitable for the two single-turbine models coupled with either the linear or the Lanzilao and Meyers merging method. Overall, the corrected Gaussian wake model with the Lanzilao and Meyers wake merging methodology gives good agreement, within 5%–6% of the LES data for all cases studied here.

Future work can be focused on studying a wider range of parameter values, such as other turbine operating conditions (larger thrust coefficients), streamwise and spanwise spacings, the ratio of the turbine hub-height to the diameter, and alternate arrangements such as partial or full staggering. This study can be extended to smooth-to-rough roughness transitions as well as other canonical configurations of surface roughness heterogeneities, such as stripes of finite widths, oblique stripes, and rectangular patches with varying aerodynamic roughness values.

The authors gratefully acknowledge the SERB Core Research Grant (CRG/2022/006735) for enabling this work as well as the computational resources made available under the National Supercomputing Mission (Grant DST/NSM/R&D_HPC_Applications/2021/28) on Param-Brahma cluster at IISER Pune and Param-Seva cluster at IIT Hyderabad.

The authors have no conflicts to disclose.

Naveen N. Kethavath: Conceptualization (equal); Data curation (lead); Investigation (lead); Writing – original draft (lead); Writing – review & editing (equal). Niranjan S. Ghaisas: Conceptualization (equal); Formal analysis (equal); Funding acquisition (lead); Writing – review & editing (lead).

The data that support the findings of this study are available from the corresponding authors upon reasonable request.

In this section, we study the sensitivity of the results to the grid size. The HT1 simulation (involving the surface roughness heterogeneity as well as the wind farm) is performed on a series of grids. The domain size is ( L x × L y × L z ) =  ( 96 d × 12 d × 4.5 d ). The aerodynamic surface roughness jumps from z 01 = 3.75 × 10 3 d to z 02 = 4.5 × 10 5 d at x s = 10 d and the first turbine row is at x t = 14 d from the left end of the main domain.

Profiles of the time-averaged streamwise velocity (U), the streamwise turbulence intensity ( T I = U U ¯ / U h), and total negative shear stress (TSS) are shown in Figs. 20–22, respectively. The results are also averaged over the three turbine columns. Vertical profiles at four downstream locations behind three rows (row 1, row 5, and row 9) are shown. Except for the coarsest mesh (case HT1-g1), the results are seen to be largely insensitive to the grid resolution with minimal differences between the two finest grids (HT1 and HT1-g5).

FIG. 20.

Mean streamwise velocity profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

FIG. 20.

Mean streamwise velocity profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

Close modal
FIG. 21.

Streamwise turbulence intensity profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

FIG. 21.

Streamwise turbulence intensity profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

Close modal
FIG. 22.

Total negative shear stress profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

FIG. 22.

Total negative shear stress profiles at different distances downstream of turbines in rows 1, 5, and 9 obtained using 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points.

Close modal

To further quantify the sensitivity to the grid, the time-disk and rotor-disk averaged streamwise velocity (U) and the streamwise turbulence intensity ( T I = U U ¯ / U h) as a function of turbine row for the five grid spacings are shown in Figs. 23(a) and 23(b), respectively. Richardson extrapolated (RE) values are determined through a mixed-order analysis,34,46 utilizing HT1-g1, HT1-g2, and HT1-g5 cases. Similar to the observations in Figs. 20–22, the profiles corresponding to the cases HT1-g2, HT1-g3, HT1, and HT1-g5 are fairly close to each other and to the profiles obtained by Richardson extrapolation. The relative errors in the mean velocity and TI, assuming the Richardson extrapolated profiles to be the exact solution, as a function of turbine row are shown in Figs. 24(a) and 24(b). For the HT1 grid spacing, the mean absolute errors are found to be 1.7% and 6.2%, respectively, for the normalized mean velocity and the turbulence intensity. Further studies for all other cases are carried out using the same resolution as in the HT1 case, i.e., with ( n x × n y × n z ) =  ( 768 × 96 × 108 ) grid points.

FIG. 23.

(a) Normalized time-averaged streamwise velocity and (b) turbulence intensity as a function of turbine row utilizing grid resolutions of 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points. RE is the result obtained from Richardson extrapolation of cases HT1-g1, HT1-g2, and HT1-g5.

FIG. 23.

(a) Normalized time-averaged streamwise velocity and (b) turbulence intensity as a function of turbine row utilizing grid resolutions of 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points. RE is the result obtained from Richardson extrapolation of cases HT1-g1, HT1-g2, and HT1-g5.

Close modal
FIG. 24.

Relative error as a function of turbine row in (a) mean streamwise velocity and (b) turbulence intensity, utilizing grid resolutions of 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points. RE is the result obtained from Richardson extrapolation of cases HT1-g1, HT1-g2, and HT1-g5.

FIG. 24.

Relative error as a function of turbine row in (a) mean streamwise velocity and (b) turbulence intensity, utilizing grid resolutions of 256 × 32 × 30 (HT1-g1), 512 × 64 × 60 (HT1-g2), 640 × 80 × 90 (HT1-g3), 768 × 96 × 108 (HT1), and 1024 × 128 × 120 (HT1-g5) grid points. RE is the result obtained from Richardson extrapolation of cases HT1-g1, HT1-g2, and HT1-g5.

Close modal

The sensitivity to the domain size is shown in Fig. 25. The HT1-g4  case employs a domain that is larger in the vertical direction. To keep the same grid spacing as in the HT1 case, the number of points, nz, is increased with respect to the HT1 case. Profiles at 3d downstream of the three turbine rows (row 1, row 5, and row 9) are seen to be largely insensitive to the vertical domain size. The vertical domain size is also in accordance with that used by previous simulations.39 

FIG. 25.

Profiles of (a) mean streamwise velocity and (b) streamwise turbulence intensity at 3d downstream of the turbines in rows 1, 5, and 9 obtained using domains of sizes ( 96 × 12 × 4.5 ) d (case HT1) and ( 96 × 12 × 6 ) d (case HT1-g4 ).

FIG. 25.

Profiles of (a) mean streamwise velocity and (b) streamwise turbulence intensity at 3d downstream of the turbines in rows 1, 5, and 9 obtained using domains of sizes ( 96 × 12 × 4.5 ) d (case HT1) and ( 96 × 12 × 6 ) d (case HT1-g4 ).

Close modal

The U and TI profiles (Figs. 26 and 27, respectively) obtained from the LES “HomogE” are validated with experimental data of Chamorro and Porté-Agel.10 Here, the boundary layer depth is δ = 4.5 d, the uniform aerodynamic roughness height is z 0 / δ = 4.444 × 10 5, and the local thrust coefficient ( C T ) values used for the different rows in the wind farm are consistent with those reported in Stevens et al.39 The LES results at different distances downstream of only four rows (row 1, row 4, row 7, and row 10) are shown for brevity, and the results are similar behind other turbine rows. The mean velocities and turbulence intensity values agree well with the experiments. Some differences in the mean velocity (Fig. 26) are seen up to 3d downstream of the first-row turbine and at 1d downstream of the other rows. The TI profiles (Fig. 27) are underpredicted behind the first-row turbine but are well predicted behind the other downstream turbine rows. These differences between the LES and the experiments may be attributed to the turbine support structure utilized in the wind tunnel measurements, which are absent from the LES. A similar level of agreement is observed in previous LES (without the turbine support structures) of this wind farm reported by Stevens et al.,39 which are also shown in Fig. 26, thus building confidence in our numerical methodology.

FIG. 26.

Mean streamwise velocities obtained from the LES case HomogE compared to the wind tunnel experiments of Chamorro and Porté-Agel10 and LES results of Stevens et al.39 Vertical profiles at four locations downstream of four representative rows are compared.

FIG. 26.

Mean streamwise velocities obtained from the LES case HomogE compared to the wind tunnel experiments of Chamorro and Porté-Agel10 and LES results of Stevens et al.39 Vertical profiles at four locations downstream of four representative rows are compared.

Close modal
FIG. 27.

Streamwise turbulence intensity profiles obtained from the LES case HomogE compared to the wind tunnel experiments of Chamorro and Porté-Agel10 and LES results of Stevens et al.39 Vertical profiles at four locations downstream of four representative rows are compared.

FIG. 27.

Streamwise turbulence intensity profiles obtained from the LES case HomogE compared to the wind tunnel experiments of Chamorro and Porté-Agel10 and LES results of Stevens et al.39 Vertical profiles at four locations downstream of four representative rows are compared.

Close modal

To ensure that the differences observed in this paper are actually caused by the changed surface roughness, rather than a different ambient turbulence intensity, we compare cases HT4, HomogT, and an additional case HomogT-z0. This last case considers a wind farm with a homogeneously rough surface with the z0 tuned such that the ambient turbulence intensity upstream of the first turbine matches that in the HT4 case. The ambient rotor-averaged TI values (i.e., 1d upstream of the first turbine row) are 0.16 for HomogT and 0.12 for HT4 and HomogT-z0. In the HomogT-z0 case, this TI is obtained by setting z 0 = 8.75 × 10 4 d. Case HT4 is chosen for this exercise rather than the cases HT1, HT2, or HT3 because the TI values 1d upstream of the first turbine differ the most between cases HomogT and HT4.

Figure 28 shows the rotor-averaged TI values 1d upstream of each turbine row, the normalized absolute and relative powers, and the wake growth rate ( k *) values at the turbine locations as a function of turbine row. The TI values upstream of the first turbine row [Fig. 28(a)] are identical between the HT4 and HomogT-z0 cases by design. However, the TI values upstream of the rest of the turbines differ slightly between the HT4 and HomogT-z0 cases. This is expected since, in the absence of the wind farm, the TI values are not expected to change with x over a homogeneously rough surface but are expected to decrease behind a rough-to-smooth transition. Figure 28(b) shows that decreasing the z0 leads to an increased absolute power for the HomogT-z0 case compared to the HomogT case. The absolute powers in the HT4 case are also very different compared to the HomogT-z0 case because of the differences in the upstream surface roughness values. The relative powers, shown in Fig. 28(c), are quite different between the HomogT-z0 and the HT4 cases. These differences are consistent with the wake growth rate values shown in Fig. 28(d). For the HomogT-z0 case, the wake recovery is slower, as characterized by the smaller k * values, which leads to larger wake losses than in the HT4 case.

FIG. 28.

(a) Turbulence intensity at 1d upstream of each turbine row, (b) normalized absolute power, (c) relative power, and (d) wake growth rate, as a function of turbine row. All values are rotor-averaged and also averaged over the three turbine columns.

FIG. 28.

(a) Turbulence intensity at 1d upstream of each turbine row, (b) normalized absolute power, (c) relative power, and (d) wake growth rate, as a function of turbine row. All values are rotor-averaged and also averaged over the three turbine columns.

Close modal

These results make it clear that the upstream ambient TI is not sufficient to characterize the wake losses in a wind farm sited behind a surface roughness heterogeneity. The streamwise variation of the background (due to the surface heterogeneity and in the absence of wind farm) velocity and turbulence intensity fields is crucial for explaining wake losses in wind farms sited behind a rough-to-smooth surface roughness heterogeneity.

1.
Abkar
,
M.
and
Porté-Agel
,
F.
, “
A new boundary condition for large-eddy simulation of boundary-layer flow over surface roughness transitions
,”
J. Turbul.
13
,
N23
(
2012
).
2.
Bastankhah
,
M.
and
Porté-Agel
,
F.
, “
A new analytical model for wind-turbine wakes
,”
Renewable Energy
70
,
116
123
(
2014
).
3.
Bastankhah
,
M.
and
Porté-Agel
,
F.
, “
Experimental and theoretical study of wind turbine wakes in yawed conditions
,”
J. Fluid Mech.
806
,
506
541
(
2016
).
4.
Bastankhah
,
M.
,
Welch
,
B. L.
,
Martínez-Tossas
,
L. A.
,
King
,
J.
, and
Fleming
,
P.
, “
Analytical solution for the cumulative wake of wind turbines in wind farms
,”
J. Fluid Mech.
911
,
A53
(
2021
).
5.
Bou-Zeid
,
E.
,
Meneveau
,
C.
, and
Parlange
,
M. B.
, “
Large-eddy simulation of neutral atmospheric boundary-layer flow over heterogeneous surfaces: Blending height and effective surface roughness
,”
Water Resour. Res.
40
,
W02505
, https://doi.org/10.1029/TR039i006p01048 (
2004
).
6.
Bou-Zeid
,
E.
,
Anderson
,
W.
,
Katul
,
G. G.
, and
Mahrt
,
L.
, “
The persistent challenge of surface heterogeneity in boundary-layer meteorology: A review
,”
Boundary-Layer Meteorol.
177
,
227
245
(
2020
).
7.
Bradley
,
E. F.
, “
A micrometeorological study of velocity profiles and surface drag in the region modified by a change in surface roughness
,”
Q. J. R. Meteorol. Soc.
94
,
361
379
(
1968
).
8.
Chamorro
,
L. P.
and
Porté-Agel
,
F.
, “
A wind-tunnel investigation of wind-turbine wakes: Boundary-layer turbulence effects
,”
Boundary-Layer Meteorol.
132
,
129
149
(
2009
).
9.
Chamorro
,
L. P.
and
Porté-Agel
,
F.
, “
Velocity and surface shear stress distributions behind a rough-to-smooth surface transition: A simple new model
,”
Boundary-Layer Meteorol.
130
,
29
41
(
2009
).
10.
Chamorro
,
L. P.
and
Porte-Agel
,
F.
, “
Turbulent flow inside and above a wind farm: A wind-tunnel study
,”
Energies
4
(
11
),
1916
1936
(
2011
).
11.
Cheng
,
W.-C.
and
Porté-Agel
,
F.
, “
A simple physically-based model for wind-turbine wake growth in a turbulent boundary layer
,”
Boundary-Layer Meteorol.
169
,
1
10
(
2018
).
12.
Du
,
B.
,
Ge
,
M.
, and
Liu
,
Y.
, “
A physical wind-turbine wake growth model under different stratified atmospheric conditions
,”
Wind Energy
25
(
10
),
1812
1836
(
2022
).
13.
Elliott
,
W. P.
, “
The growth of the atmospheric internal boundary layer
,”
Eos, Trans. Am. Geophys. Union
39
,
1048
1054
(
1958
).
14.
Ghaisas
,
N. S.
, “
A predictive analytical model for surface shear stresses and velocity profiles behind a surface roughness jump
,”
Boundary-Layer Meteorol.
176
,
349
368
(
2020
).
15.
Ghaisas
,
N. S.
,
Ghate
,
A.
, and
Lele
,
S. K.
, “
Effect of tip spacing, thrust coefficient and turbine spacing in multi-rotor wind turbines and farms
,”
Wind Energy Sci.
5
,
51
72
(
2020
).
16.
Ghate
,
A. S.
and
Lele
,
S. K.
, “
Subfilter-scale enrichment of planetary boundary layer large eddy simulations using discrete Fourier-Gabor modes
,”
J. Fluid Mech.
819
,
494
539
(
2017
).
17.
Ghate
,
A. S.
,
Ghaisas
,
N. S.
,
Lele
,
S. K.
, and
Towne
,
A.
, “
Interaction of small scale homogeneous isotropic turbulence with an actuator disk
,” AIAA Paper No. 2018-0753,
2018
.
18.
Goit
,
J. P.
and
Önder
,
A.
, “
The effect of coastal terrain on nearshore offshore wind farms: A large-eddy simulation study
,”
J. Renewable Sustainable Energy
14
,
043304
(
2022
).
19.
Howland
,
M. F.
,
González
,
C. M.
,
Martínez
,
J. J. P.
,
Quesada
,
J. B.
,
Larranaga
,
F. P.
,
Yadav
,
N. K.
,
Chawla
,
J. S.
, and
Dabiri
,
J. O.
, “
Influence of atmospheric conditions on the power production of utility-scale wind turbines in yaw misalignment
,”
J. Renewable Sustainable Energy
12
(
6
),
063307
(
2020
).
20.
Howland
,
M. F.
,
Lele
,
S. K.
, and
Dabiri
,
J. O.
, “
Wind farm power optimization through wake steering
,”
Proc. Natl. Acad. Sci. U. S. A.
116
(
29
),
14495
14500
(
2019
).
21.
Ishihara
,
T.
and
Qian
,
G. W.
, “
A new Gaussian-based analytical wake model for wind turbines considering ambient turbulence intensities and thrust coefficient effects
,”
J. Wind Eng. Ind. Aerodyn.
177
,
275
292
(
2018
).
22.
Jimenez
,
A.
,
Crespo
,
A.
,
Migoya
,
E.
, and
García
,
J.
, “
Advances in large-eddy simulation of a wind turbine wake
,”
J. Phys.: Conf. Ser.
75
,
012041
(
2007
).
23.
Katic
,
I.
,
Højstrup
,
J.
, and
Jensen
,
N. O.
, “
A simple model for cluster efficiency
,” in
European Wind Energy Association Conference and Exhibition
(
A. Raguzzi
,
Rome
,
Italy
,
1986
), Vol.
1
, pp.
407
410
.
24.
Kethavath
,
N. N.
,
Mondal
,
K.
, and
Ghaisas
,
N. S.
, “
Large-eddy simulation and analytical modeling study of the wake of a wind turbine behind an abrupt rough-to-smooth surface roughness transition
,”
Phys. Fluids
34
(
12
),
125117
(
2022
).
25.
Lanzilao
,
L.
and
Meyers
,
J.
, “
A new wake-merging method for wind-farm power prediction in the presence of heterogeneous background velocity fields
,”
Wind Energy
25
(
2
),
237
259
(
2022
).
26.
Lissaman
,
P. B. S.
, “
Energy effectiveness of arbitrary arrays of wind turbines
,”
J. Energy
3
(
6
),
323
328
(
1979
).
27.
Mondal
,
K.
,
Kethavath
,
N. N.
,
Abhinay
,
K.
, and
Ghaisas
,
N. S.
, “
Large eddy simulation study of atmospheric boundary layer flow over an abrupt rough-to-smooth surface roughness transition
,”
Boundary-Layer Meteorol.
188
,
229
257
(
2023
).
28.
Mulhearn
,
P. J.
, “
Relations between surface fluxes and mean profiles of velocity, temperature and concentration, downwind of a change in surface roughness
,”
Q. J. R. Meteorol. Soc.
103
,
785
802
(
1977
).
29.
Niayifar
,
A.
and
Porté-Agel
,
F.
, “
Analytical modeling of wind farms: A new approach for power prediction
,”
Energies
9
(
9
),
741
(
2016
).
30.
Panofsky
,
H. A.
and
Townsend
,
A. A.
, “
Change of terrain roughness and the wind profile
,”
Q. J. R. Meteorol. Soc.
90
,
147
155
(
1964
).
31.
Peña
,
A.
,
Schaldemose Hansen
,
K.
,
Ott
,
S.
, and
van der Laan
,
M. P.
, “
On wake modeling, wind-farm gradients, and AEP predictions at the Anholt wind farm
,”
Wind Energy Sci.
3
(
1
),
191
202
(
2018
).
32.
Porté-Agel
,
F.
,
Bastankhah
,
M.
, and
Shamsoddin
,
S.
, “
Wind turbine and wind farm flows: A review
,”
Boundary-Layer Meteorol.
174
,
1
59
(
2020
).
33.
Rao
,
K. S.
,
Wyngaard
,
J. C.
, and
Coté
,
O. R.
, “
The structure of the two dimensional internal boundary layer over a sudden change of surface roughness
,”
J. Atmos. Sci.
31
,
738
746
(
1974
).
34.
Roy
,
C. J.
, “
Grid convergence error analysis for mixed-order numerical schemes
,”
AIAA J.
41
(
4
),
595
604
(
2003
).
35.
Rozema
,
W.
,
Bae
,
H. J.
,
Moin
,
P.
, and
Verstappen
,
R.
, “
Minimum-dissipation models for large-eddy simulation
,”
Phys. Fluids
27
,
085107
(
2015
).
36.
Schulz-Stellenfleth
,
J.
,
Emeis
,
S.
,
Dörenkämper
,
M.
,
Bange
,
J.
,
Cañadillas
,
B.
,
Neumann
,
T.
,
Schneemann
,
J.
,
Weber
,
I.
,
Zum Berge
,
K.
,
Platis
,
A.
et al, “
Coastal impacts on offshore wind farms–A review focussing on the German Bight area
,”
Meteorol. Z.
31
,
289
315
(
2022
).
37.
Shir
,
C. C.
, “
A numerical computation of air flow over a sudden change of surface roughness
,”
J. Atmos. Sci.
29
,
304
310
(
1972
).
38.
Stevens
,
R. J. A. M.
and
Meneveau
,
C.
, “
Flow structure and turbulence in wind farms
,”
Annu. Rev. Fluid Mech.
49
,
311
339
(
2017
).
39.
Stevens
,
R. J.
,
Martínez-Tossas
,
L. A.
, and
Charles
,
M.
, “
Comparison of wind farm large eddy simulations using actuator disk and actuator line models with wind tunnel experiments
,”
Renewable Energy
116
,
470
478
(
2018
).
40.
Stevens
,
R. J.
,
Gayme
,
D. F.
, and
Meneveau
,
C.
, “
Coupled wake boundary layer model of wind-farms
,”
J. Renewable Sustainable Energy
7
(
2
),
023115
(
2015
).
41.
Subramaniam
,
A.
,
Ghate
,
A.
,
Ghaisas
,
N. S.
,
Howland
,
M. F.
et al, PadeOps GitHub Repository, see https://github.com/FPAL-Stanford-University/PadOps/tree/igridSGS (last accessed August 5, 2021).
42.
Teng
,
J.
and
Markfort
,
C. D.
, “
A calibration procedure for an analytical wake model using wind farm operational data
,”
Energies
13
(
14
),
3537
(
2020
).
43.
Troen
,
I.
and
Lundtang Petersen
,
E.
, European Wind Atlas,
1989
, see https://backend.orbit.dtu.dk/ws/portalfiles/portal/112135732/European_Wind_Atlas.pdf.
44.
Vahidi
,
D.
and
Porté-Agel
,
F.
, “
A physics-based model for wind turbine wake expansion in the atmospheric boundary layer
,”
J. Fluid Mech.
943
,
A49
(
2022
).
45.
Van Der Laan
,
M. P.
,
Pena
,
A.
,
Volker
,
P.
,
Hansen
,
K. S.
,
Sørensen
,
N. N.
,
Ott
,
S.
, and
Hasager
,
C. B.
, “
Challenges in simulating coastal effects on an offshore wind farm
,”
J. Phys.: Conf. Ser.
854
,
012046
(
2017
).
46.
van der Laan
,
M. P.
and
Abkar
,
M.
, “
Improved energy production multi-rotor wind farms
,”
J. Phys.: Conf. Ser.
1256
,
012011
(
2019
).
47.
Voutsinas
,
S.
,
Rados
,
K.
, and
Zervos
,
A.
, “
On the analysis of wake effects in wind parks
,”
Wind Eng.
14
(
4
),
204
219
(
1990
).
48.
Wu
,
Y. T.
and
Porté-Agel
,
F.
, “
Large-eddy simulation of wind-turbine wakes: Evaluation of turbine parametrisations
,”
Boundary-Layer Meteorol.
138
(
3
),
345
366
(
2011
).