Evaluation of Wind Flow Characteristics by RANS-Based Numerical Site Calibration (NSC) Method with Met-Tower Measurements and Its Application to a Complex Terrain

The performance of wind turbines is not only dependent on the wind turbine design itself, but is also dependent on the accurate assessment of wind resources at the installation site. In this paper, the numerical site calibration (NSC) method using three-dimensional Reynolds-averaged Navier–Stokes (RANS) simulation was proposed to accurately forecast the wind flow characteristics of wind turbine sites with complex terrains, namely Methil in Scotland, and Haenam in South Korea. From NSC at the Methil and Haenam sites, it was shown that the complicated and vortical flow fields around hills and valleys were captured using the three-dimensional RANS CFD simulation in Ansys CFX software based on a high-resolution scheme with a renormalization group (RNG)-based k-ε turbulence model. It was also shown that topographically induced wind profile and turbulence intensity over a local-scale complex terrain are remarkably dominated by flow separation after passing hills. It was concluded that the proposed NSC method using three-dimensional RANS simulation with a high-resolution scheme was an economically useful method for evaluating wind flow characteristics numerically to assess wind turbine sites with complex terrains and designing the wind farm layout.


Introduction
Wind energy has been in the spotlight as a major source of renewable energy, and its levelized cost of energy (LCOE) has reduced a lot to become a favorable alternative source to fossil fuel energy. With innovatively engineered design of blade airfoils and planform using increased computational performance, wind turbines have been designed to extract wind energy as close as possible to the Betz limit of 59.3% [1].
Since onshore wind turbines are installed not only in flat areas as preferred sites, but also in areas surrounded by hills or valleys, it is important to understand the geological effect on wind flow approaching the wind turbine to maximize the annual energy production (AEP). Therefore, in order to accurately evaluate the performance of a wind turbine overall, it is essential to forecast the wind flow characteristics such as velocity profile and turbulence intensity related to the wind turbine installation site. Furthermore, it is recommended to choose a wind turbine of a design grade adequate for the estimated wind conditions in order to maximize durability and economic benefits [2][3][4].
There are many studies investigating the effect of wake on wind speed and the performance of a single wind turbine or multiple wind turbines like on a wind farm. Ainslie outlined the features of wind turbine wakes and described a numerical model to calculate the wake flow field [5]. Vermeer et al. studied the aerodynamics of horizontal axis wind turbine wake and reviewed both the near and the far wake regions. For the far wake, both a single turbine and wind farm effects were experimentally and numerically studied to understand how the far wake decays downstream [6]. Sanderse et al. also reviewed numerical calculations of wind turbine wake aerodynamics with different computational fluid dynamics (CFD) techniques for modeling the rotor and the wake [7]. Several works were focused on the validation of wake simulation models with measurements from wind tunnel experiments and directly from wind farms. Chamorro and Porte-Agel performed wind tunnel experiments to study turbulence in the wake of a wind turbine model [8], and Schumann et al. investigated wind turbine wakes in a wind tunnel with two similar models of wind turbine [9]. Gaumond et al. used power production data from the Horns Rev and Lillgrund offshore wind farms to compare three engineering wake models [10]. In addition, Barthelmie et al. evaluated the commonly used wake models for predicting wind speed decrease downstream of a wind turbine based on six experiments using a ship-mounted sodar at a small wind farm [11]. Specifically, Barthelmie et al. used integrated Doppler lidar data from the Perdigao experiment to measure wind turbine wakes and also to characterize wind gusts in a complex terrain [12,13]. Energy extraction out of various wind farm layouts was investigated with a wake effect model [14], and it showed at most a 30% decrease in the extracted energy depending on the wind farm layouts. Jadhav and Roy [15] examined the impact of wind turbine wake on the economic dispatch of a power system and stated that the wake effect significantly influenced the overall cost of operation and the loss of transmission. Wind turbine loads under the wake effect were studied based on the measurements obtained from Vindeby offshore wind farm in Denmark [16], which showed that the turbulence intensity and fatigue loading in the wind turbine can increase by 5-15% according to the layout of wind farm. Kim et al. [17] also investigated the wake effect from the wind farm layout and showed a 7% reduction in AEP and a 30-50% increment in damage equivalent load (DEL).
In addition, there are many interesting works related to noise caused by a complex terrain. A major obstacle to the diffusion and worldwide spread of wind energy is represented by its produced noise [18,19]. Wind turbine do not emit high levels of noise, but they are much more annoying than other common noise sources [20,21]. Especially in hilly areas, quietness is required and even lower noises are perceived much more. Prevention and monitoring play key roles in avoiding well-known health related issues such as sleep disorders [22], learning impairment [23,24], hypertension ischemic heart disease [25], and annoyance [26]. In this regard, a proper knowledge of winds will also be useful to better predict the noise emission of a wind farm through detailed procedures [27].
According to the International Electrotechnical Commission (IEC) 61400-12-1 [28], site calibration procedures must be followed before installing wind turbines except for regions that show only minor variations in topography. For a site calibration, anemometers are installed at the prospective wind turbine site and a reference site nearby, and the correlation of the wind conditions at the two sites would be derived based on the data acquired during a specific period. The correlation is then used to forecast the wind characteristics at the prospective installation site. For this correlation process, two meteorological towers have to be erected to measure the speed and the direction of the wind, requiring significant amounts of time and cost. Moreover, the site calibration is suitable for a relatively flat topography with stable winds, and therefore it is difficult to obtain accurate forecasts in topographically complicated regions such as Korea, which has a lot of mountainous terrains resulting in unstable winds.
According to Annex B of the IEC 61400-12 standard, the site calibration is required except for cases where the maximum slope value and the maximum terrain variation from the plane at the prospective installation site are smaller than the recommended specific values. The standard suggests that a flow model can be used to determine whether a site calibration measurement can be avoided or not if the terrain characteristics are within an additional 50% of the limits of maximum slope value. In that case, if the difference in wind speed between the anemometer and the installation site is within 1% at 10 m/s, then the site calibration measurement is not necessary.
Sanada et al. [29] conducted studies on forecasting unstable winds at complex terrains using a numerical site calibration method (NSC) based on flow models from computational fluid dynamics (CFD) [30,31]. The cross-validation of the numerical simulations with a high-resolution scheme and k-ε model against wind tunnel experimental results and on-site measurements indicates a good agreement [32].
In this paper, the numerical site calibration method (NSC) employing three-dimensional Reynolds-averaged Navier-Stokes (RANS) approach was proposed to numerically forecast the wind flow characteristics at the wind turbine installation sites without the installation of expensive meteorological towers for the type certification of a single wind turbine. Two specific areas, Methil in Scotland and Haenam located to the west of South Korea, were first investigated with the given topography data and a three-dimensional RANS CFD simulation solver in Ansys CFX [33].

Numerical Site Calibration (NSC) Method
The numerical site calibration method (NSC) is a three-dimensional flow analysis method based on computational fluid dynamics (CFD) to forecast the wind condition, and it utilizes the correlation of wind characteristics between the reference sites and the prospective wind turbine installation sites. Compared to traditional site calibration, NSC eliminates the cost of erecting a meteorological tower and can be applied to complex terrains, providing wind resource to accurately estimate the performance of large wind turbines as shown in Figure 1 [34]. In order to establish the evaluation method based on NSC, it is important to clarify the range of application for calculation results and its reliability. Other points which should be also addressed include the characteristics of the numerical analysis code, calculation regions and boundary conditions, and securing statistical significance. Since calculation regions and boundary conditions are important factors for simulating wind conditions, wind condition comparison between measurement data and simulation values must be performed. The comparison needs to show the correlation error between measured and simulated values within 3% to establish the simulation procedure.

Governing Equation and Turbulence Flow Model
There are three major numerical techniques for modeling turbulent flow fields: direct numerical simulation (DNS), large eddy simulation (LES), and Reynolds-averaged Navier-Stokes (RANS). DNS and LES require substantial calculation resources to analyze a vortex, and it is difficult to apply them to practical engineering problems. With recent developments in computer capacity and capability, there have been a number of studies on the specific behavior of a vortex by performing numerical analyses that apply LES on complex terrains [35][36][37]. On the other hand, RANS uses time-based, ensemble-averaged Navier-Stokes equations and models all of the effects from turbulence. Although RANS yields a lower resolution of analysis than DNS or LES, it is widely used in engineering applications due to the practical aspect of not requiring high-resolution calculation grids [38][39][40][41].
Unlike the complex flow fields inside and outside a wind turbine, in most cases of wind condition analysis, it is difficult to form rapid adverse pressure gradients in the flow direction. Therefore, the k-ε turbulence model, which accurately analyzes turbulence behavior in free-stream regions with small pressure gradients, is chosen preferably for RANS numerical analysis than the Wilcox k-ω turbulence model [42], which provides accurate forecasts of boundary layer separation having the adverse pressure gradients in viscous sub-layer regions. The k-ε turbulence model is widely used by researchers for topography analysis [43]. In particular, the renormalization group (RNG)-based k-ε turbulence model yields forecasts with high reliability in analyzing complex terrains [44]. For this reason, the RNG-based k-ε turbulence model was used for this paper. Even though wall grid-scale is not enough near the surface, the current numerical analysis was also tested using the shear stress transport (SST) turbulence model, which is based on the k-ω turbulence model in order to evaluate the sensitivity of the turbulence model [45]. In this paper, the k-ε turbulence model has been used though two approaches resulted in little difference. The convergence of numerical analysis was determined based on the regular wind speed history at the meteorological tower site.
In the RANS equation, the flow is considered as the ideal gas equation of state. For 3D steady flow simulation, the instantaneous equations of mass, momentum and energy conservation in tensor form can be written as follows [45]. The continuity, momentum and energy equations are given by In Equations (1)- (3), ρ is the air density, u j is the jth velocity component, τ ij is the Reynolds stress tensor of eddy viscosity model, p is the air static pressure, e is the air energy, C p is the specific heat, µ is the air viscosity, Pr is the Prandtl number, and T is the air temperature. The air viscosity and Prandtl number are given by where µ l is the laminar viscosity, µ t is the turbulent viscosity calculated by a turbulence model, Pr l is the laminar Prandtl number, and Pr t is the turbulent Prandtl number. Version 15.0 of Ansys CFX was used in this simulation.

Wind Turbine Sites
In this paper, two specific sites, Methil in Scotland, and Haenam in South Korea, were selected to evaluate the wind flow characteristics using the proposed NCS method based on RANS simulation as shown in Figures 2 and 3. In both cases, WAsP software was used to transform the Google map into digitalized information for CFD simulation. The Methil site is a sort of flat area with a few steep hills and residential areas, while the Haenam site is quite a complex terrain consisting of many hills and valleys, which is a common characteristic in South Korea.  The velocity at the Met tower at a height of 80 m was measured for six months and Table 1 shows 10 min average velocities at specific wind directions, where the main wind direction was 210 deg.

Boundary Conditions
To perform numerical site calibration properly, calculation regions and boundary conditions should be defined. Figure 4 shows the calculation regions for this analysis which were set at 5.0 km with a minimum grid interval of 5 m in all directions. In order to minimize the blockage effect on the analytical results, a sufficient distance of 5 km was maintained from the ground in the vertical direction. As per the boundary conditions, a no-slip condition was applied on the ground surface in the vertical direction, and the roughness length of the sea, land, and residential area were set at 0.0002, 0.03, and 0.4 m, respectively, as written in IEC 61400 part 12 ANNEX A [28]. In addition, a free-slip condition was applied to the ceiling of the calculation region. The wind shear profile of Equation (6) was applied as both entry and exit boundary conditions in the east-west-north-south directions.
In Equation (6), V and H denote wind speed and vertical height, respectively, and the subscript hub means the position of the wind turbine hub. α is a heuristic constant that varies according to the atmospheric conditions. In this study, α = 0.2 was used based on the IEC standard. In this work, a turbine with 171 m rotor diameter and 110 m hub height was used. Table 2 lists the summary of the boundary conditions for numerical site calibration at the Methil site.
Here, k and ε are computed by the turbulence transportation equation, and the outlet condition was computed through the pressure boundary. Since the Mach number is very low, incompressibility can be assumed, therefore temperature fluctuation is ignorable. The reason for using the RNG k-ε model is that it is difficult to simulate the first cell with a very small height, but it has more freedom in selecting the y + parameter. The simulation convergence was good enough since rms pressure was converged below 10 −6 .

Grid Sensitivity Study
It is necessary to understand the characteristics of the flow analysis code to evaluate the reliability of the calculation outcome. In this study, the effects of grid densities in the horizontal and vertical directions and the turbulence model on the calculation results were evaluated. Table 3 summarizes the number of cells throughout the calculation region. The grid density distributions in the horizontal direction of grid system numbers 1 and 2 in Table 3 are shown in Figure 5, which indicates that grids were formed around the wind turbine and the meteorological tower. For the current mesh generation, the fully hexahedral mesh was used.  Table 3. Grid systems and analysis cases.

Methil, Scotland
The digitalized raw image was pre-processed by ICEM to generate grid meshes, then analyzed by Ansys CFX. The numerical results were post-processed using Fieldview software for the high-resolution image as shown in Figure 6. As expected, the overall topography shows little variation in terms of elevation with only scattered hills. Based on the topography data, the maximum slope value and the maximum terrain variation from the plane were calculated according to the IEC standard part 12 [28]. The calculation results exceeded the specified values due to the high hill on the planned location of the meteorological tower, then it was required to follow the site calibration procedure by installing an additional Met tower at the wind turbine installation location. The numerical analysis based on the 3D RANS simulation was performed to find out the resolution to avoid the required site calibration process.
From the numerical analysis, the effect of the height of the hill on the wind speed difference between the tower location (reference site) and the prospective wind turbine installation site was investigated. The difference in wind speed between the reference site and the prospective wind turbine installation site in the identical wind direction of 210 deg at a boundary wind speed condition of 10.2062 m/s is summarized in Table 4. As listed in Table 4, the effect of grid density on the calculation result is just 0.13%, which is ignorable. Furthermore, in order to quantitatively evaluate the effect of the turbulence model on the calculation result during RANS simulation, sensitivity analysis was performed on the k-ε turbulence model and the shear stress transport (SST) turbulence model, which has the characteristics of the k-ω turbulence model. Table 5 indicates that the effect of the turbulence model on calculation result is very minimal at 0.05%. Therefore, it could be stated that the uncertainties associated with the effects of grid density and the turbulence model on the calculation outcome are very little.    Table 6 shows results of the difference in wind speed between the two sites according to the decrease in the height of the hill where the meteorological tower is to be installed. If the height of the hill decreases to 6 m, the difference in wind speed converges to within 1%, which can potentially eliminate the need for site calibration measurement. Table 7 and Figure 8 also show that the difference in wind speed is within 1% across the valid measurement sectors if the hill height is adjusted to 6 m, which confirms that site calibration measurement could be eliminated.   Figure 9 shows velocity error based on the varying Met tower terrain altitudes and the close-up images of the area within a 753 m radius (4.4 × rotor diameter) of the wind turbine, which confirms that as the height of the hill decreases from the original location (21.95 m) to 6 m, the difference in wind speed between the reference site and the prospective installation site decreases within the acceptable range of 1% described in IEC standard.

Haenam, South Korea
3D RANS simulation was performed on the region of Songho-ri, Haenam-gun, Jeollanam-do, Republic of Korea. It is a coastal area with complex mountainous terrains. Identical procedures were followed as for Methil, Scotland for establishing calculation regions and performing analyses. Figure 10 displays the elevation distribution of the topography in Haenam, which indicates a narrow strip of land and a maximum elevation of 400 m above sea level. Compared to Methil, where the maximum elevation is 50 m with a wider land area, Haenam has much more variation in terms of elevation. Such a complex terrain causes a large-scale vortex, which increases turbulence intensity and the fatigue load on the wind turbine. Therefore, 3D RANS CFD analysis was performed to assess the potential risks of the prospective wind turbine installation site by forecasting the turbulence intensity.
From the numerical results as shown in Figure 11, it visualizes the wind speed distribution of the equal-interval cross-section of a complex terrain and areas with turbulence intensity of 15% or greater (indicated in black). Because the k-ε turbulence model was applied to the 3D RANS simulation in this study, turbulence intensity (T i ) can be calculated using Equation (7).
In Equation (8), u and U denote the varying velocity component and velocity, respectively. In Equation (9), k represents turbulence kinetic energy. As shown in Figure 11, in the southern mountainous terrain, the wind speed is high at regions with high turbulence intensity. In the northeastern region, however, the wind speed is low where the turbulence intensity is high. The discrepancy was caused by the 3D vortex structure formed in the complex terrain, which showed that there was little correlation between wind speed and turbulence intensity in a complex terrain. Figure 12 displays a close-up image of the vortex flow field of the northeastern region (box area in white dot) of the Haenam site shown in Figure 11. Figure 12a displays the flow speed distribution in the equal-interval cross-section, regions with high turbulence intensity, and the flow around the region with high turbulence intensity. The analytical results of Figure 12a indicate that Region 1 (white dotted line) has high wind velocity and low turbulence intensity, whereas Region 2 has low wind speed and high turbulence intensity. Unlike Region 1, flows of vortex structure have formed in areas around Region 2. Figure 12b displays the elevation distribution in a complex terrain, limiting streamline on the ground level, and areas with high turbulence intensity. A tornado-shaped vortex structure is formed in the area indicated with the white dotted line in Figure 12b, which indicates high turbulence intensity of more than 15%. Therefore, if a wind turbine is installed near Region 1 of Figure 12a or in an area with high turbulence intensity in Figure 12b, there can be performance degradation and durability problems caused by fatigue load [46].

Discussion
The numerical site calibration method based on 3D RANS simulation was proposed to forecast and evaluate wind flow characteristic speeds, and this was demonstrated at two specific sites, a simple terrain with steep hills (Methil in Scotland) and a complex terrain with hills and valleys (Haenam in Korea).
The current method to forecast the wind flow characteristics is an applicable and economically useful way to eliminate the additional Met tower installation. In the case of Methil, by adjusting the hill height to 6 m, the estimated value is very close to the measured wind data. However, artificially adjusting hill height from 21.95 m to 6 m is quite a difficult task including resident permission, environmental damage, etc. For a flat region like Methil, however, the estimated wind speed is still reasonable with 3.75% deviation from the measured value without adjusting hill height. Therefore, the application of the current method requires close site investigation including topological conditions.
In the Haenam site with a complex terrain, high turbulence intensity of more than 15% can be caused by the 3D vortex structure, which could result in performance degradation and durability problems due to fatigue load. In this case, special attention should be given to a wind turbine system to relieve fatigue loads or safety execution to protect the system from unwanted fatigue loading.
The results of the current study can be summarized as follows.

1.
Steeper hills in Methil significantly affect the wind flow at the measurement height of the wind turbine and the meteorological tower. Based on CFD analysis results, if the installation altitude of the meteorological tower is lower (6 m) than the original altitude (21.5 m), costly site calibration measurements can be avoided.

2.
The topography-induced wind profile and the turbulence intensity over local-scale complex terrains such as valleys and hills are dominated by separation flows in Haenam. With the proposed numerical method, the risks to turbine performance and durability could be significantly avoided by assessing the wind characteristics of the wind turbine installation site.

3.
The proposed numerical site calibration method based on three-dimensional RANS simulation is very useful for evaluating wind flows at prospective wind turbine installation sites over complex terrains and could be utilized for optimizing the layout of a wind turbine array.