**1. Introduction**

Limited carbon resources and environmental concerns are some of the reasons leading the energy industry to exploit alternative energy sources. Wind energy systems have been developed and applied for sites with suitable conditions, the first modern commercial-scale wind turbines were placed in United States approximately 40 years ago. Nowadays, the most common and profitable applications for wind energy systems are the large wind farms. Commercial-scale wind generators for wind farms are within 3 MW and 5 MW, and all have a predominantly horizontal axis and are three bladed. One problem of these large wind farms is the row arrangemen<sup>t</sup> of the generators. The towers are usually placed in rows, requiring large areas of land for rotors up to 100 m in diameter. Previous research has suggested safe distances to avoid the wind turbines blade/components damage and output power waste. However, the optimum spacing between turbines in a wind farm is still a challenging and open question in wind energy research.

Several efforts using different methodologies have been done to achieve layout optimization, focusing on finding optimal spacing between turbines in a wind farm. Park and Law [1] applied sequential convex programming to maximize wind farm output power by optimizing the placement of wind turbines of the Horns wind farm in Denmark. They found that the optimal spacing between wind turbines is dependent on the wind direction. Scattering the turbines helped to avoid wake chain effects, so that downstream rotors were not significantly affected. Moreover, the same study considered

wind statistical data to optimize the wind farm power production over a long period, resulting in a 7.3% power increase. Son et al. [2] found that the total wind farm output power is strongly related to the distance between the first and second wind turbine rows. When the referred distance became larger, the output power considerably dropped in comparison to smaller distances. This means that the increase of the spacing between the first and second rows is ineffective in improving output power. On the other hand, decreased distances made the second wind turbine row much less efficient. They discovered the importance of keeping turbines as close as possible, but with enough space so that the second row can have guaranteed output power. Longer distances did not contribute to increase the total output power. Further, increasing the space between the fourth and the fifth rows has a better contribution than increasing the space between the first and second rows. Wu and Porté-Agel [3] investigated two layout configurations in the same area with 30 turbines either arranged in aligned or staggered conditions. In comparison to the aligned configuration, the staggered one allows better wake recovery. This exposes the downstream turbines to higher local wind speeds (consequently higher performance) and lower turbulence intensity. Stevens et al. [4] found that the distance of 10 diameters (or higher) would minimize the cost per unit of energy production, and the same is true for a distance of 15 diameters if the objective function was evaluated using dimensionless parameters. Those value are significantly higher than applied values in wind farms (6–10 turbine diameters). Meyers and Meneveau [5] found that the current wind farms layout solutions in literature have characteristics with considerably lower spacing than computationally optimized layout solutions.

Moreover, other efforts have attempted to achieve wind farm optimization using control strategies to mitigate wake effects, applying sub-optimal operating conditions. This means that each rotor will not necessarily deliver the best aerodynamic performance, but the goal is to find the best solution that avoids wake interaction effects, increasing the total wind farm output power. Park and Law et al. [6] studied control strategies for wake effects mitigation, showing that control techniques can be applied for each individual rotor to improve overall wind farm efficiency. González et al. [7] proposed the individual selection of an operating point on each wind turbine in order to maximize the overall wind farm output power. This is performed by studying the optimal pitch angle and Tip Speed Ratio (TSR) of each rotor in regards to the total wind farm output power. Additionally, the methodology also allows decreased turbulence intensity levels in the produced wakes. The results showed increased power production when the wind speed is lower than the rated wind speed, and for non-prevailing wind directions. Lee et al. [8] found an increase of 4.5% in the total output power by applying pitch angle control for the Horns Rev wind farm. Kazda et al. [9] applied weakened wake conditions for upstream turbines by using sub-optimal operations through control strategies. They found that a 12.5% reduction for the upstream turbines resulted in a 2.5% increase in the sum of the upstream and downstream turbines. This could be achieved by either a change of 3.5◦ in the pitch angle or by a 24% reduction in TSR compared to optimum TSR. For the case of two upstream turbines operating at 87.5% of optimal conditions, the sum of total power of the upstream and downstream turbines increased by 9.7%. Gil et al. [10] applied control strategies, achieving from 1.86% up to 6.24% in energy captured by using sub-optimal operating points. Chowdhurry et al. [11] found that using variable rotor diameters improved efficiency, achieving 30% increase in the total power generation.

All these efforts in the literature described above provided relevant contributions to wind farm optimization and turbine spacing research. However, they did not consider a rigorous evaluation of three-dimensional wake effects, which this study will achieve. Most of the computational studies from literature proposed design optimization frameworks in which the wind turbine models were not based on data validated against experimental measurements of the wake flow field. In the context of science applied to wind farm optimization, this work proposes a numerical Computational Fluid Dynamics (CFD) model to rigorously analyze wind turbine wake flow field, characterizing wake flow characteristics for the most relevant wind farm parameters: velocity flow field and turbulence intensity. The current study will do a full computational analysis of the near-wake aerodynamic behavior, considering configurations not analyzed before in literature: several different loading, free-stream velocity and pitch angle conditions. The goal is to achieve a validated model by comparing computational and experimental data from existing literature. Engineering tools such as CFD or wake analytical methods have been improved to accurately characterize wake characteristics, but there are few experiments to effectively validate wind turbine wake flow. The present study provides an advance in relation to previous models as it propose a physical model in which the near wake is validated against a critical wind farm design parameter: the velocity flow field. With a reliable CFD model validated in terms of flow field, this model could be applied to improve wind farm layout design by including both near and far wake regions. Literature shows a variety of techniques and different goals in regards to wind turbine CFD. The next section shows a description of the main experimental approaches found in literature, which will be useful to provide data to develop and validate the numerical model in this study.

#### *1.1. Brief Description of Wind Tunnel Experiments*

A full review of low-speed wind tunnel studies and scaled turbines is provided by Crespo et al. [12]; additionally, other recent relevant studies can be found in the literature [13–22]. Most of these studies are meant to validate wind turbine simulations, and some of them are described below to provide an overview of low-speed wind tunnel experiments. The objective of this literature review is to show the way that experimental data can be used in order to validate wind turbine simulations.

Wind turbine experiments conducted by the Norwegian University of Science and Technology validated the numerical results against wind tunnel measurements in terms of mean velocity, turbulence intensity and the power and thrust coefficients. This research center has low-speed wind tunnel facilities, with dimensions of 2.71 m wide, 1.8 m high and 11.1 m long. An experimental study was performed using two aligned prototype rotors of 0.944 m and 0.894 m, and the blade consists of 14% S826 NREL (National Renewable Energy Laboratory) profile for the two rotors [13]. The velocity profile was characterized using Pitot-Static tubes, and the thrust force was determined using a Six-Component Balance Force.

A qualitative study of the rotor wake behavior by Chamorro and Porté-Angel [14] analyzed a 150 mm diameter three-bladed wind turbine prototype, which was tested using a wind tunnel with 37.5 m length driven by a 200 hp fan. The experimental data was used to produce a qualitative study of the wake behavior, since the Reynolds number is different compared to full-scale wind turbines. A particularly interesting aspect that distinguishes this study from the others is that the authors were able to characterize the surface roughness by placing straight chains of approximately 5 mm height covering a 10 m section of the tunnel. These chains were aligned perpendicular to the flow direction and separated from each other by 0.20 m. The mean wind velocity in the tunnel was measured using Pitot static tubes, and constant TSR values (λ = 4.2 for smooth surfaces and λ = 4.4 rough surfaces) were maintained in order to reflect the typical operational conditions of full-scale field turbines (typically 3 < λ < 6).

In another experiment, a virtual wind-tunnel model (24.4 m × 36.6 m) with the same dimension of the NASA (National Aeronautics and Space Administration) wind tunnel was analyzed using the ANSYS Fluent package [15]. The model validation was performed comparing the pressure coefficient at different span-wise sections along the turbine blade. In addition, the wind turbine output power was compared to published experimental results for the NREL phase VI rotor tested in the NASA wind tunnel. Several other studies in literature utilized data from the NREL/NASA framework to develop CFD studies using pressure coefficient values on the blades and aerodynamic torque data for comparison and validation. Zhou et al. [23] performed LES (Large Eddy Simulation) of the NREL phase IV to evaluate the effect of different inflow conditions on the aerodynamic loading and near wake characteristics. Hsu et al. [24] implemented a finite-element (Lagrangian-Eulerian) model of the NREL Phase IV using a non-structured rotating mesh refined close to the rotor disc. Wake characterization was not the focus of the study, what explains the wake made out of coarse non-structured cells with no refinement. Gundling et al. [25] evaluated low and high fidelity models using the NREL Phase

VI for predicting wind turbine performance, aeroelastic behavior and wakes: 1) The Blade Element method with a free-vortex wake; 2) The actuator disc method; 3) The full-rotor method. Mo et al. [26] did a study in more depth to understand wake aerodynamics performing a LES of the NREL Phase VI using dynamic Smagorinsky-model, additionally verification of the average Turbulence Intensity was performed against an analytical model. They found that the downstream distance where instability and vortex breakdowns occur is dependent on wind free-stream inlet conditions (7 <sup>m</sup>·s<sup>−</sup><sup>1</sup> happens at four rotor diameters, while 15.1 <sup>m</sup>·s<sup>−</sup><sup>1</sup> between 11 and 13 diameters), and a decrease of the turbulence intensity happened after instability and vortex breakdowns. Choudhry et al. [27] performed a very similar CFD study of the NREL phase VI using computational methods similar to the ones found in the study conducted by Mo et al. [26], finding that regions of velocity deficit and high turbulence intensity are within the high vorticity region.

Sturge et al. [16] utilized an open-circuit suction tunnel, driven by an eight-blade axial fan positioned at the outlet. In this experiment, the wind speed is controlled by using a variable frequency drive. The air flow passes through a honeycomb mesh with cells 0.01 m wide and 0.1 m long. The dimensions vary along the tunnel, with a 6.25:1 contraction section and 1.2 m high × 1.2 m wide × 3 m long test section. Afterwards, analysis of static pressure along the blade showed a large reduction in the suction peak along the leading edge, which reduced the lift generated by the rotor and consequently the torque production.

The wake flow of a 5 × 5 array of 50 mm micro-wind turbines was studied and analyzed by Houssain et al. [18] using a wind tunnel. These 1/10 scaled prototypes were placed ina3m × 1.8 m wind tunnel, allowing the velocity profile and turbulence intensity (velocity fluctuations) behind the array to be measured at different downstream locations. The wake flow was characterized by using hot-wire anemometer, ultrasonic anemometer measurements, and Particle Image Velocimetry (PIV). The full-scale rotor of 500 mm diameter was analyzed as well. The results for velocity deficit and the turbulence intensity were similar for both rotors.

In this sense, the MEXICO (Model Experiments in Controlled Conditions) experiment [28] was one of the most comprehensive collaborative efforts by the International Energy Agency (IEA), who created the task 29 to gain understanding about wind turbine aerodynamics, as well as to improve aerodynamic models used for wind turbine design. A series of tests for a small wind turbine prototype were performed using the DNW German Dutch open section wind tunnel. Although the rotor wake measurements comprised only the near wake region right behind to the wind turbine (up to 1.33D downstream of the rotor), the experiment is a very rich source of data useful to validate wind turbine CFD wake models.

This present work covers the gap of characterizing the wind turbine wake flow field based on experimental data from existing literature, which describes the validation of a wind turbine CFD simulation using velocity wake data from the MEXICO experiment. The goal is to extend the understanding of the wake flow field beyond the distances analyzed in these experiments, and also analyzing the influence of variable operating conditions on near wake aerodynamic behavior. In order to do so, variable operating conditions with regards to the TSR and the Pitch Angle (θ) were simulated to understand how these specific design parameters affect the flow field. The second part of this work will extend the analysis beyond the near wake, characterizing the far wake aerodynamic behavior according to the same TSR and Pitch Angle (θ) conditions.

#### *1.2. Detailed Overview of the MEXICO Experiment*

The experiments described in the previous section only performed rotor measurements. However, computational models based on CFD assumptions also need flow field measurements to be successfully validated. The most comprehensive experimental flow field measurement study was the MEXICO Experiment [28], which used a rotor prototype of 4.5 m diameter and the largest wind tunnel existent in the European continent. PIV techniques were employed to collect flow field measurements around

the rotor plane (Figure 1). Several recent studies utilized data from the MEXICO experiment to validate their CFD models [29–45] with different research goals as detailed below.

**Figure 1.** Sketch showing an overview of the MEXICO (Model Experiments in Controlled Conditions) Experiment (Top View).

In regards to Lifting Line codes, Yang et al. [29] showed the necessity for developing new techniques to account for 3D rotational effects on predicting loading for rotors. They created a new technique to determine the angle of attack on rotating blades using data from the MEXICO experiment, a Blade Element Momentum (BEM) code relying on 2D airfoil data was found to over-predict the loading of the rotor; this discrepancy was attributed to the 3D effects originated from the rotor geometry. Xudong et al. [30] developed an aerodynamic/aero-elastic design tool to optimize wind turbine blades and validated the results using MEXICO data for turbine loading.

Regarding the first round of PIV wake measurements (axial flow), Bechmann et al. [31] performed a CFD simulation of the MEXICO rotor using RANS (Reynolds Averaged Navier Stokes) equations further downstream up to 2.5 diameters behind the rotor. All the simulations were done fully turbulent, but there might be laminar flow at the leading edge of the blades; further work is needed to demonstrate the length of accuracy of laminar turbulent-transition models. Micallef et al. [32] characterized the radial velocities in the near wake close to the MEXICO rotor using a potential-flow panel model to characterize the wake radial induction. Tip vortex characterization performed by tracking its location showed that the radial flow velocity in the rotor plane is not fully dominated by the blade vorticity. Carrión et al. [33] assumed periodic boundary conditions to model only one of the MEXICO rotor blades under axial flow conditions, finding good agreemen<sup>t</sup> for the wake flow field by using a compressible multi-block solver without needing to switch between compressible and incompressible flow. Herraez et al. [34] validated a CFD model in OpenFoam using the Spalart-Allmaras turbulence model, showing comparisons for pressure distributions from several blade sections, and PIV near wake measurements. Shen et al. [35] performed CFD simulations of the MEXICO rotor including the geometry of the wind tunnel, and regarding tunnel wall effects this study found that tunnel effects are not significantly influenced by the fluid flow. Garcia et al. [36] developed a hybrid filament-mesh vortex method to improve computational efficiency, using the MEXICO experimental dataset for near wake validation. Nilsson et al. [37] described vortex structures in the near wake of the MEXICO rotor using the actuator line method. The trajectory of the tip vortices and wake expansion were described according to the TSR, implementing a RANS LES model. Wimshurst and Willden [38] simulated the near wake flow field of the MEXICO rotor using multiple reference frame approach. The actuator line method using 2D aerodynamic data was compared to a 3D polar actuator line model. Zhong et al. [39] developed a numerical tool combining Lagrangian dynamic large-eddy and actuator line models using PIV wake data for validation, finding that the tip vortices contribute to a maximum velocity deficit peak and turbulence intensity peak near the blade tip. Guntur and Sørensen [40] developed a full rotor CFD model of the MEXICO rotor focusing on the flow at the inboard part of the blades, analyzing the

boundary layer separation at this region to understand differences in behavior between 3D flow and 2D flow. This latter study showed that the fluid flow separation starts at a higher angle of attack for the 3D case.

In regards to the second round of measurements (yawed flow), Sørensen et al. [41] did the first attempts to validate the near wake flow field in yawed flow. Tsalicoglou et al. [42] performed RANS computations of the MEXICO rotor wake for yawed and uniform flow cases, showing that the velocity deficit in the near wake (up to two diameters downstream) does not follow a Gaussian distribution. Additionally, the interaction with structures of the wind turbine (nacelle and tower) is more significant for yawed flows. The effects on the wake caused by the tower and the blade could still be observed at the end of the near wake. Grasso and Garrel [43] showed that the lifting line code coupled with the free wake method can accurately represent the near wake at uniform or yawed conditions. Shen et al. [44] developed an actuator line/Navier-Stokes model using the MEXICO rotor experimental dataset under yawed flow for flow field validation, considering both loading and velocity flow field for the simulation.

## **2. Computational Methods**

#### *2.1. Rotor Blade Geometry*

The MEXICO experiment performed several different flow field measurements to characterize the three-dimensional velocity flow field in the near wake. Experimental measurements such as traverse and longitudinal wind velocity, both upwind and downwind of the rotor, were performed at a few specific locations. Here, we validated the computational model by plotting the velocity in the wake region of the blade and directly comparing the simulation results with experimental data from the MEXICO rotor. Because our hope is to implement a rapid computational simulation, the objective is to obtain agreemen<sup>t</sup> between experimental and computational velocities within 5%. The rotor simulated in this work was the MEXICO Rotor (Figure 2); the three-bladed model has three types of airfoil: DU91-W2-250 (20% to 45%), Riso-A1-21 (54% to 65%), and NACA 64-418 (75% until the blade tip). The blade is also twisted, and a pitch angle of −2.3◦ was applied for the measurements. The blade geometry can be found in the final report of this experiment [28]. Since some of the airfoil data are not publicly available, a reverse engineering process was performed to find the airfoil coordinates.

**Figure 2.** MEXICO rotor geometry, a three-bladed rotor with 4.5 m diameter. Reference for the blade geometry data: Scheppers et al. [28].

#### *2.2. Layout and Boundary Conditions*

We broke down the computational domain (Figure 3) into smaller parts for two reasons. First, local mesh sizing: the meaningful region can be refined to correctly characterize the flow field. Second, pressure-far-field boundary conditions for the lateral and superior boundaries require a larger domain to keep straight streamlines at the boundaries to achieve numerical convergence. The dimensions of the square part containing the wind tunnel and the rotor extends from −2.5D to 2.5D, while the exterior part corresponding to the surroundings extends from −10D to 10D.

**Figure 3.** Layout of the computational domain and boundary conditions. (**a**) Perspective view and boundary conditions; (**b**) Information about the physical domain; (**c**) Front view of the central disc.

#### *2.3. Computational Fluid Dynamics Modeling*

CFD assumptions are based on the Finite Volume Method (FVM) for representing and evaluating partial differential equations in the form of algebraic equations. The domain of interest is divided into small cells, reducing the Navier-Stokes equations to algebraic or simple differential equations. Integration of the volume is conducted to obtain surface fluxes because the flux entering a given volume is identical to that leaving the adjacent volume. The CFD solver implemented in this work was ANSYS Fluent 17, housed in two computers, each with 64 GB RAM/8 processes. ANSYS Fluent solves the equations of fluid flow and heat transfer by default using a stationary (or inertial) reference frame. However, a Moving (or non-inertial) Reference Frame (MRF) can bring advantages in solving the equations for some problems involving moving parts, such as rotating blades. In those problems, the flow around the moving parts is the variable of interest to be determined. In the case of this work, the region behind of the wind turbine corresponding to the wake flow field is the region of interest. The MRF technique models the flow around the moving part as a steady-stead problem with respect to the moving frame, allowing to activate reference frames in selected cell zones. The ANSYS Fluent MRF modeling modify the equations of motion to incorporate additional acceleration terms that occur due to the transformation from the stationary to the moving reference. The main reason for employing a MRF is to solve a problem that is unsteady in the stationary (inertial) frame but steady with respect to the moving frame. In this work, the simulation was performed using a steady state MRFapproach, and setting the rotational speed to match experimental conditions. The turbulence model selected was the k–ω SST, which is suitable for swirl flow, and it was used in the literature

studies as their main turbulence modelling technique. Because there is no public information from the reports of the MEXICO experiment regarding the inlet inflow conditions, default values of 5% and 10% were assumed for the inflow turbulence intensity (TI) and the viscosity ratio (VR) at the inlet, respectively. The Reynolds number based on the average chord length is approximately 1.5 × 105. Pressure-far-field boundaries are applied for the lateral and superior boundaries, pressure-outlet for the exit, velocity-inlet for the front boundary, and a special type of wall with no shear (named No Shear Wall inferior) for the inferior boundary, which represents the bottom of the physical domain (Figure 3). Different operating conditions were tested in this experiment, and some of them were mimicked in this computational study for the validation: ω = 424.5 rpm, U = 15 <sup>m</sup>·s<sup>−</sup><sup>1</sup> (which results in a TSR = λ = 6.6), and U = 10 <sup>m</sup>·s<sup>−</sup><sup>1</sup> (TSR = 10). Additionally, several other operating conditions regarding Free-Stream Velocity, TSR and Pitch Angle were simulated to characterize the wake aerodynamic behavior.

The physical domain was meshed using unstructured elements (Figure 4), which are suitable for CFD applications because of its good convergence rate. The mesh sensitivity study showed a total of approximately 10 million cell elements to be sufficient to accurately validate the model and describe the near wake (Appendix A). The meshing process consisted of a sphere of influence with 0.1m cell elements in a radial distance of 6 m surrounding the rotor, and a square part extending from −0.5D to 3D with 0.25 m cell elements. The blade surface mesh was dimensioned using local edge sizing to reduce the skewness of the cells, resulting in 175 nodes spanwise and 75 nodes chordwise at the blade tip. Additionally, 10 inflation layers with a ratio of 1.1 were built to ensure y+ < 1 next to the blade surface. The physical domain needs to be large enough to result in a good simulation convergence, since pressure-far-field boundaries (lateral boundaries) require straight streamlines to avoid divergence for the residuals. However, the mesh at the exterior part surrounding the wind turbine and the rotor domain is coarse, since this region is not meaningful for the CFD analysis.

**Figure 4.** *Cont.*

**Figure 4.** Mesh of the Computational Domain. (**a**) Computational domain; (**b**) Rotational disc; (**c**) Front view: sectional plane showing details of the rotational central disc and the sphere of influence; (**d**) Lateral view: sectional plane showing details of the sphere of influence and the rectangular near wake region; (**e**) Sectional plane showing the mesh close to the blade surface; (**f**) Sectional plane showing inflation layers close to the blade surface.

#### *2.4. TSR (λ) Effect on the Near Wake*

A very important design parameter for wind farms is the TSR, which is defined as the ratio between the blade tip speed velocity and the free-stream velocity (Equation (1)). The TSR and other parameters such as free-stream velocity are critical to determine the wake behavior:

$$
\lambda = \frac{\boldsymbol{\omega} \cdot \mathbb{R}}{\mathbb{U}\_{\text{freestream}}} \tag{1}
$$

where ω is the rotor rotational speed, R is the blade radius and U is the free stream velocity.

Another important design parameter is the Turbulence Intensity (TI). This parameter can be calculated using the Equation (2):

$$\text{TI} = \frac{\sigma\_{\text{U}}}{\text{U}\_{\text{freestream}}} \tag{2}$$

where σU is the velocity standard deviation.

#### *2.5. Wake Validation*

The flow field at the wake of the rotor is validated by comparison between experimental [28] and computational data from the CFD simulation. The axial and radial traverses at the wake described in the section 1.2 (Figure 1) are considered for the validation. The MEXICO experimental dataset is an extensive one, with different turbine configurations (axial, yawed, azimuth angle) under multiple operating conditions (velocity, TSR, stand-still). In this work, specific conditions were selected for the validation: axial flow for U = 10 <sup>m</sup>·s<sup>−</sup><sup>1</sup> and U = 15 <sup>m</sup>·s<sup>−</sup>1. This decision is related to the complexity of the dataset, which would make it challenging to mimic using a steady-state CFD model. Rigorously, a CFD model would not only need to replicate velocity flow field but also tip vortex tracking agreement. There is no such a work in literature, the vast majority of the studies selected specific operational conditions and turbine configurations in order to narrow the scope of the research. Even though a wind speed of 24 <sup>m</sup>·s<sup>−</sup><sup>1</sup> is available in the MEXICO experimental dataset, different near wall resolutions may have a considerable effect on the prediction the wake flow field. Particularly, the turbulence profile in the wake could be significantly affected by different orders of near wall resolution. As pointed out by Shen et al. [35], higher free stream wind speeds would require a more refined mesh for a more accurate prediction of the wake flow field.

#### *2.6. Near Wake Analysis*

Besides implementing and validating the CFD model of the MEXICO rotor, simulations were carried out considering variable operating conditions other than the ones analyzed in the original experiment. As a result of such simulations, a detailed study on near wake aerodynamics was developed to assess numerical sensitivity of the wake in regards to TSR, Pitch Angle and Free-Stream Velocity.
