A multiscale numerical framework coupled with control strategies for simulating a wind farm in complex terrain

A multiscale numerical framework coupled with control strategies for simulating a wind farm in complex terrain

Energy 203 (2020) 117913 Contents lists available at ScienceDirect Energy journal homepage: www.elsevier.com/locate/energy A multiscale numerical f...

3MB Sizes 0 Downloads 0 Views

Energy 203 (2020) 117913

Contents lists available at ScienceDirect

Energy journal homepage: www.elsevier.com/locate/energy

A multiscale numerical framework coupled with control strategies for simulating a wind farm in complex terrain Qiang Wang , Kun Luo *, Renyu Yuan , Shuai Wang , Jianren Fan , Kefa Cen State Key Laboratory of Clean Energy Utilization, Zhejiang University, Hangzhou, 310027, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 13 July 2019 Received in revised form 25 April 2020 Accepted 17 May 2020 Available online 20 May 2020

Improving the accuracy of the evaluation on the performance of wind farms in complex terrain under the actual atmosphere is crucial to the sustainable development of wind power. To this end, a multiscale numerical framework (MNF) integrating the WRF model coupled with a wind farm parameterization and a LES model embedded with the wind turbine control strategies is developed. An interface for data exchange between the WRF and LES models is established. The simulated nacelle velocity and the generator power results agree well with the observations, proving that the MNF is capable of assessing the flow behavior and performance of a realistic onshore wind farm in complex terrain. Besides, the impacts of terrain on the operating process for wind turbines over the hill under the dynamic inflow condition are quantified. The results show that the wind turbine at the top of the hill can cut in the ratedstate earlier and cut out later, which improves the power output by 5.5%, indicating that the total generation of wind farms can be enhanced by making full use of terrain acceleration. This high-fidelity multiscale numerical framework is expected to be an effective tool for the siting and optimization of wind farms. © 2020 Published by Elsevier Ltd.

Keywords: Wind farm Complex terrain Multiscale numerical framework Control strategy Operating process

1. Introduction With the rapid development of the wind power industry, a large number of flat terrain wind resources have been utilized. Nowadays, wind developers pay more attention to exploit the wind resources in complex terrain [1]. For instance, the distributed wind application in a complex environment has become one of the most essential policies in China [2]. However, assessing the operation characteristics of wind farms accurately, in particular for those in the real atmosphere, is still challenging. Therefore, there is an urgent need to develop a high-fidelity numerical framework for simulations of the wind farm in complex terrain. In general, the approach for assessing the wind farm involves in situ observation, tunnel test, and numerical simulation [3]. The in situ observation is a direct way to investigate the dynamic wake flows and the operating processes of the wind farm in natural wind conditions by using measuring equipment. It includes the traditional anemometry devices [4] and some innovative devices, e.g., Light Detection and Ranging (LiDAR) based on an alternative

* Corresponding author. E-mail address: [email protected] (K. Luo). https://doi.org/10.1016/j.energy.2020.117913 0360-5442/© 2020 Published by Elsevier Ltd.

mobile measurement technology [5] and the instrumented drone [6]. However, the in situ observation mainly focuses on the discrete locations, and the uncertainty arises from the irregular variability, randomness, and intermittent of the atmosphere due to the limitation of the observation technique. To reduce the uncertainty of the test, the scaling wind turbines were made to test the wake properties using the wind tunnel. However, the most tunnel experiments carried out in uniform or empirical inflow boundary condition for understanding the steady flow situation (see the review by Chi-Jeng Bai [7]). In summary, though they could serve as the basic data for numerical model validation, it is not advisable to study the wind farm only use experimental methods given the respective limitations, e.g., expensive and time consuming. As an alternative approach, the numerical method can flexibly obtain comprehensive and continuous information of the wind farm in complex terrain and atmospheric boundary layer (ABL) conditions. Besides, it is not restricted by the spatial and temporal scales of the wind farm as well. According to the scale of the wind farm, the numerical simulations can be divided into two categories: (i) the mesoscale simulation using the numerical weather prediction (NWP) models and (ii) the microscale simulation using the computational fluid dynamics (CFD) models. The mesoscale simulation is used to assess the wind farm based

2

Q. Wang et al. / Energy 203 (2020) 117913

Nomenclature

Variables A An, Bn cv CP CT CTKE D eL eD kn K KI KP IDr Ls n NGear P t T Taero TGen U Ua Uhub Ut v10,ave z

rotor area Fourier coefficient vectors inflow speed correction constant power coefficient thrust coefficient turbulence kinetic energy coefficient rotor diameter unit vectors of lift force unit vectors of drag force wavenumber vector roughness coefficient of the ground pitch controller integral gain pitch controller proportional gain drivetrain inertia of the shaft turbulence integral scale unit sphere vector high-speed to low-speed gearbox ratio power simulation time inflow speed fluctuation period aerodynamic torque generator torque wind speed axial inflow speed speed at hub-height relative tangential speed average speed at a height of 10 m vertical height

Greek symbols angle of attack random unit sphere vector

a an, bn

on the relatively coarse grid resolution by using a numerical weather prediction model. Reviewing the literature, it can be seen that the majority of mesoscale simulations focused on the model sensitivity analysis, the wind farm parameterization method development, and the applications of wind resource and power prediction of wind farms. Siuta et al. [8] and Mughal et al. [9] systematically evaluated the sensitivity of hub-height wind speed to the planetary boundary layer (PBL) scheme, grid length, and initial condition selection in the WRF model over complex terrain, indicating that the maximal improvements in simulation accuracy primarily depended on the grid length or PBL scheme. Dzebre et al. [10] and Carvalho et al. [11,12] tested the sensitivity of surface winds to 11 of the PBL schemes available in the WRF model, showing that the MYNN scheme can accurately reproduce the local wind resource characteristics. Furthermore, to represent the wind farm in the realistic ABL, the wind farm parameterization (WFP) concept has been proposed by Baidya Roy [13], treating the wind farm as a momentum sink and turbulence source. Then, the WFP was successively developed and optimized by Blahak et al. [14], -Agel et al. [16]. Recently, a list of relevant Fitch et al. [15], and Porte nez et al. [17] reproduced the power studies has been reported. Jime deficits associated with wind turbine wakes in an offshore wind farm. Lee et al. [18,19] investigated the wind farm wake evolution during the evening transition, showing that the presence of the wind farm would increase the wind speed deficit monotonically

l h hGB ufn un zf r q qk 4

U Urate

persistence parameter Kolmogorov scale gearbox efficiency natural frequency Fourier frequency damping ratio density of air blade pitch angle corrected blade pitch angle inflow angle rotor rotation speed rated rotor rotation speed

Abbreviations ABL Atmosphere Boundary Layer ALM Actuator Line Method BEM Blade Element Momentum CFD Computational Fluid Dynamics DWT Downwind Wind Turbine GK Gain-correction Factor GTC Generator-Torque Control KS Kinematic Simulation LES Large-Eddy Simulation LiDAR Light Detection and Ranging MNF Multiscale Numerical Framework NWP Numerical Weather Prediction PIC Proportional-Integral Control PBL planetary boundary layer TP Time length Proportion TSR Tip Speed Ratio TWT Top Wind Turbine UDF User-Defined Functions UWT Upwind Wind Turbine WFP Wind Farm Parameterization WPS WRF Preprocessing System WRF Weather Research and Forecasting

and make surface cooling. Mangara et al. [20] used the WFP scheme to simulate turbulent wake produced by a real onshore wind farm. In summary, the mesoscale WRF model can gain the wind characterisitcs and performance of the wind farm in the real ABL, but it is difficult to accurately represent the flow details and the operating process for each wind turbine in wind farm. On the other hand, the microscale simulation could resolve the blade-element or turbine scale turbulence by using the CFD models, which plays a crucial role in accurately modeling wind farms, including wake dynamics, operating characteristics, wind farm siting, and power prediction. Currently, there are three CFD modeling approaches [21], i.e., direct numerical simulation (DNS), Reynolds-averaged Navier-Stokes (RANS), and large eddy simulation (LES). Recently, LES is being widely used by institutions and enterprises for scientific research and business practice and has become a high-fidelity means for the development and utilization -Agel of wind farms (see the review by Mehta et al. [22] and Porte et al. [23]). In microscale simulations, it is reasonable to model each wind turbine as an actuator in the wind farms [24,25] to reduce computational cost. However, the inflow boundary condition is a challenge to microscale simulations. Current wind farm simulations are mostly initialized by using the empirical logarithmic or power inflow wind profile [26,27]. It can be concluded that the microscale simulations, in particular for LES, can gain high-fidelity simulations for wind farms, but the experienced inflow boundary

Q. Wang et al. / Energy 203 (2020) 117913

conditions need to be improved to be more realistic. Detecting the real-like ABL condition is an indispensable assignment to improve wind farm modeling accuracy. Considering the merits of each model, the combination of mesoscale and microscale models has become a promising approach for evaluating the performance of wind farms. Among them, the most common way is to conduct a CFD simulation initialized by using the results obtained from the mesoscale model.  et al. [28] analyzed the wake behaviors of two aligned wind Mache turbines using a coupled WRF-LES model, showing that the turbulence becomes more homogeneous after the second turbine. Carvalho et al. [29] applied the WRF model with the WAsP software to assess the wind of a region in complex terrain and the output of an onshore wind farm. Gopalan et al. [30] developed a coupled mesoscale and microscale model to estimate the wind resource and aerodynamics of wind turbines with good capability. However, the previous combination of mesoscale and microscale models still has some limitations: i) the inflow boundary conditions generated by the results from mesoscale simulations did not take the impact of large-scale wind farms into account, and ii) the optimal control strategy of wind turbine under dynamic inflow conditions has not been considered in the microscale model. Therefore, there is still room for improving the coupling framework and the accuracy for simulating the wind farms in complex terrain. Therefore, the objective of this study is to develop a multiscale numerical framework (MNF) integrating the WRF model coupled with a wind farm parameterization and the LES model embedded with the wind turbine control strategy. An interface for data exchange between the methods in different scales is constructed. After validation, this multiscale numerical framework is applied to investigate the flow behavior and performance of a realistic onshore wind farm in complex terrain and quantify the impact of

3

terrain on the operating property under dynamic inflow condition. The rest is structured as follows: the methodologies of the multiscale numerical framework are elaborated in Section 2. The description of a realistic onshore wind farm in complex terrain located in Inner Mongolia of China and the model configurations is given in Section 3. Then, the results of model validation, flow behavior and performance, and operating process characteristics of the wind turbines in wavy terrain under variable inflow condition are quantitatively analyzed in Section 4. In the end, the conclusions are presented in Section 5.

2. Methodology 2.1. Multiscale numerical framework In the current multiscale numerical framework for simulating realistic wind farms, the WRF model coupled with a wind farm parameterization is combined with the LES model embedded with wind turbine control strategy through an interface for data exchange between the mesoscale and microscale models. To bridge the gap, the mesoscale data is firstly fitted to generate a mean inflow profile. Then, turbulence fluctuation based on the local wind resource and terrain is generated and imposed by using the Kinematic Simulation (KS) method [31,32]. Consequently, a real-like inflow condition can be generated for the LES model. This implementation process is outlined in Fig. 1, including the following three steps. Step 1: Preprocessor of the wind profile in mesoscale. The regional wind field is simulated by the mesoscale WRF model based on the ARW solver coupled with the WFP model, which is

Fig. 1. Schematic diagram of the multiscale numerical framework.

4

Q. Wang et al. / Energy 203 (2020) 117913

driven by the WRF preprocessing system (WPS) and initialized by the assimilation data (e.g., FNL). Step 2: Interface of data exchange between the mesoscale and microscale. Fitting the prepared coarse wind data properly as a mean inflow wind profile, then imposing the turbulence fluctuation based on the local wind resource and the terrain on this mean inflow wind profile based on the KS method. As a result, the fine wind profile is prepared for the microscale simulations. Step 3: Running the main model for simulating the target wind farm. In the main model, the LES solver is driven by the meshed geometric model, i.e., the wind turbine modeled as ALM, and the physical model representing the turbine operating properties, in particular for the coupled control strategies like generatortorque control (GTC) and proportional-integral control (PIC) [33]. The flow behavior and operating process characteristics of the target wind farm can then be evaluated by examining the variables from the data post-processor. It should be noted that a platform of open-source CFD toolbox OpenFOAM [26] is used for microscale simulations in the present study. The remainder of this section will introduce the major pieces in Fig. 1. 2.2. WRF model coupled with WFP To prepare the coarse wind field toward the target wind farm in the innermost of the regional area, a preprocessor of the wind profile in mesoscale is firstly provided. The WRF model coupled with a wind farm parameterization (WFP) is used to reproduce the interaction between the wind farm and ABL (Step 1 in Fig. 1). Here, the WRF model and the coupled WFP are described as below. The Advanced Research WRF (ARW) dynamic solver is a stateof-the-art atmospheric prediction system in the WRF model, it integrates fully compressible, non-hydrostatic Euler equations, which are formulated using a terrain-following mass vertical coordinate. The scientific and algorithmic approaches in the ARW solver, including the governing equations, model discretization, turbulent mixing and model filters, initial conditions, lateral boundary conditions, nesting, physics, and data assimilation (see the NCAR technical note [34]). The ARW solver has the ability to resolve the flow around the wind farm in a real atmospheric environment, which is used to obtain an initial wind condition for the target wind farm simulation. A wind farm parameterization (WFP) developed by Fitch et al. [35] is used to reproduce the interaction between the wind farm and the ABL. The WFP scheme represents the impacts of wind turbines on the local atmosphere by adding a momentum sink and a turbulence source. It assumes that the total energy absorbed by a wind turbine from the atmosphere could be partly converted into useful electrical power, and partly into the turbulence kinetic energy. For brevity, the core equations and the relative explanations for the WFP model are listed in Table 1. For details please refer to the previous study [15]. The momentum tendency term in Eq. (1) can be represented by the thrust coefficient CT. The electrical power in Eq. (2) and the turbulence kinetic energy in Eq. (3) with the area of Aijk can be resolved by the power coefficient CP and turbulence kinetic energy coefficient CTKE, respectively. Where the |U|ijk is the velocity magnitude in the grid cell (i, j, k), and zkþ1 and zk represent the lower and upper-level heights of the grid cell with vertical index k.

Table 1 The core equations for the wind farm parameterization.

vjUjijk vt

¼

2 1 CT jUjijk Aijk 2 ðzkþ1  zk Þ

(1)

3 vPijk 1 CP jUjijk Aijk ¼ 2 ðzkþ1  zk Þ vt

(2)

3 vTKEijk 1 CTKE jUjijk Aijk ¼ 2 ðzkþ1  zk Þ vt

(3)

and microscale should be prepared for the microscale simulation (Step 2 in Fig. 1). Firstly, fitting the prepared coarse wind data properly as a mean inflow wind profile, and then imposing the turbulence fluctuation based on the local wind resource and the terrain on this mean inflow wind profile based on the Kinematic Simulation method [32]. The KS method generates the fluctuation speed at the inlet based on the classical Davenport spectrum, which agrees well with the wind speed field in the actual atmospheric boundary layer. In the KS model, the turbulence is simulated by a large number of random Fourier modes varying in space and time. The fluctuation speed at the inlet is calculated as follows:

u0 ðx; tÞ ¼

N X ½An cosðkn , x þ un tÞ þ Bn sinðkn , x þ un tÞ

(4)

n¼1

where N is the number of Fourier modes, An and Bn are the Fourier coefficient vectors, kn is the wavenumber vector, and un is the Fourier frequency. The wavenumber vector kn is given by:

kn ¼ kn n; kn ¼ k1

 ðn1Þ=ðN1Þ Ls

h

(5)

where n is a random unit sphere vector-centered at the point (0 0 0), Ls is the turbulence integral scale in the flow field, and h is the Kolmogorov scale. In the wavenumber spectrum, the lowest wavenumber and the highest wavenumber are defined as:

k1 ¼

2p ; Ls

kN ¼

2p

(6)

h

The selection of the Fourier coefficient vectors An and Bn needs to satisfy a principle that is:

An , kn ¼ Bn ,kn ¼ 0

(7)

where An and Bn can be achieved by:

A n ¼ An

kn  an ; jkn  an j

Bn ¼ Bn

kn  bn jkn  bn j

(8)

In this expression, an and bn are the random unit sphere vector. The Eqs. (7) and (8) make sure that the fluctuation speed resulting from Eq. (4) satisfies the continuity equation [36]. The magnitudes of the Fourier coefficient vectors An and Bn are set as:

2.3. Kinematic Simulation method

A2n ¼ B2n ¼ 2Eðkn ÞDkn

After obtaining the coarse wind data from the mesoscale simulation, an interface of data exchange between the mesoscale

where E (kn) is the value of the energy spectrum, and Dkn can be calculated as:

(9)

Q. Wang et al. / Energy 203 (2020) 117913

8 k2  k1 > > n ¼ 1; > > 2 > > > < Dkn ¼ knþ1  kn1 n ¼ 2; …; N  1; > 2 > > > > > k  kN1 > N : n¼N 2

the actuator elements.

hε ¼ (10)

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi (11)

where l is the persistence parameter [37]. For the Davenport spectrum, the value of the energy spectrum is calculated as [38]:

Eðkn Þ ¼

4Kv210;ave kn

  2  r exp  ε ε3 p3=2 1

,h

ð1200kn =2pÞ2 1 þ ð1200kn =2pÞ

i 2 4=3

(15)

where ε is the projection width, and r is the distance from the actuator point to the point where the force is applied. Meanwhile, Troldborg et al. [44] suggest that ε ¼ 2Dx, where Dx is the grid size near the actuator point. Thus, the body force for each actuator point of the rotor is expressed as:

In addition, the Fourier frequency un in Eq. (4) is given by:

un ¼ l k3n Eðkn Þ

5

(12)

where K is the roughness coefficient of the ground, v10,ave is the average speed at a height of 10 m above ground. Accordingly, the wind speed fluctuation generated by the current method depends on the terrain roughness and the vertical distribution of wind speed.

f ðx; y; z; tÞ ¼ f 5 hε ¼

N X

f ðxi ; yi ; zi ; tÞ

i¼1

  2  r exp  i ε ε3 p3=2 1

(16) where (xi, yi, zi) is the coordinate of the ith actuator point, and N is the total number of actuator points. To accurately capture the interaction between the wind turbine and the surrounding flow field, the lift and drag forces can be determined by several basic parameters of the selected wind turbine, including the airfoil design parameters (i.e., airfoils type, chord, twist, and thickness), tower, and nacelle data, which can be obtained from the wind power manufacturer. Besides, the wind turbine control strategies, such as GTC and PIC, embedded in the microscale model are introduced as below.

2.4. LES model coupled with the control strategy After preparing the fine wind profile based on the coarse data from the mesoscale simulation, the LES is performed and used to resolve the wind farm (Step 3 in Fig. 1). Here, the methodologies of the LES model with ALM, control strategies such as GTC and PIC are briefly described as below. 2.4.1. LES with actuator line model The large-eddy simulation (LES) solver is used to resolve the filtered incompressible Navier-Stocks equations by means of a pseudo-spectral discretization in the horizontal directions and a second-order finite central differencing scheme in the vertical direction [25,39]. The sub-grid scale stress is closed by the Smagorinsky model. The wind turbine induced body force is modeled using the ALM method [40e42], which is proposed by Sørensen and Shen [43]. In this model, the wind turbine blades are represented by a set of actuator points along each blade axis based on the BEM theory. The blade forces can be determined by the design blade parameters including airfoil type, chord, twist, and lift/drag coefficients. The relative speed is a dominant factor determining the magnitude of forces, which is given by:

Urel ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi U 2a þ ðUr þ Ut Þ2

(13)

where U represents the rotor speed, Ua represents the axial inflow speed, and Ut represents the relative tangential speed. The vector relationship of the above speeds is determined by the blade pitch angle q, the inflow angle 4, and the angle of attack a. The body force of each actuator element produced by the wind flow can be calculated by:

1 2 f ¼ ðL; DÞ ¼ rUrel cðCl eL þ Cd eD Þ 2

(14)

where eL and eD respectively stand for the unit vectors of lift and drag force. While the effect of flow fields on the rotor is reflected by mapping the force relations. Here, the Gaussian function is utilized to represent the body force projection on the flow field surrounding

2.4.2. Control strategy model In order to keep the maximum energy conversion efficiency, the wind turbine has to operate at a constant tip speed ratio [45]. Therefore, to improve the wind turbine operating efficiency and reduce the load on a wind turbine, the operating state should be timely adjusted according to the local inflow wind conditions. When the inlet wind speed is smaller than the rated wind speed, it is necessary to adjust the generator torque so that the wind turbine can keep a constant tip speed ratio. If the inlet wind speed is larger than the rated speed, it is vital to vary the pitch angle to avoid the wind turbine rotor speed exceeding its rated value. This is the generator-torque control (GTC) strategy. In this controller, the rotor speed U is primarily dominated by the rotor aerodynamic torque and the generator torque, and the rate of rotor speed change can be obtained by:

dU TAero  hGB  TGen  NGear ¼ dt IDr

(17)

where TAero is the rotor aerodynamic torque, hGB is the gearbox efficient, TGen is the generator torque, NGear is the gearbox ratio between the high-speed shaft and the low-speed shaft, and IDr is the drivetrain inertia of the low-speed shaft. The aerodynamic torque is driven by the rotor, and the generator torque is transmitted from the generator via the high-speed shaft. The generator torque is a piecewise function of the generator speed n, which can be determined by the generator speed in five regions [46], termed as Region 1, Region 1-1/2, Region 2, Region 2e1/2, and Region 3, as shown in Fig. 2. The turbine will operate based on the optimal torque-speed curve when the rotor is working at an optimal TSR. However, the actual torque-speed curves vary as the generator speed increases. The generator torque in Region 1 is zero when speed is less than the cut-in value. In Region 1-1/2, the generator torque undergoes a linear transition star-up stage. Then, the generator torque rises along the optimal curve, which is proportional to the square of the generator speed to keep an optimal TSR. To limit the TSR at rated value for reducing the noise of wind turbines, another linear transition Region 2e1/2 should be

6

Q. Wang et al. / Energy 203 (2020) 117913

Fig. 2. Generator torque-speed response of the variable-speed controller used in the present study.

experienced [33]. In Region 3, as the speed exceeds the rated value, the power output keeps constant, and as a result, the generator torque is inversely proportional to the generator speed. Since the generator has a minimum start-up speed and the generator torque cannot increase indefinitely, the actual torque-speed curve does not agree with the optimal curve. To control the rotor and generator speed at the rated condition, the blade pitch needs to be adjusted in real-time. Here, the PIC is introduced. It assumes that the adjustment signal of the blade pitch angle is treated as a linear combination of the proportional signal and integral signal for speed deviation. Thus, the blade pitch angle signal Dq can be obtained by:

(18)

0

where KP and KI are the proportional and integral gains for the PIC, respectively, and DU is the rotor speed deviation relative to the rated speed. Here, the gains of KP and KI are given by:

8 2IDr U0 zf ufn > > > KP ¼ N < Gear ðvP=vqÞ > > > :K ¼ I

IDr U0 u2fn NGear ðvP=vqÞ

(19)

where U0 stands for the rated rotor speed, zf stands for the damping ratio, and ufn stands for the natural frequency. Moreover, vP/vq represents the sensitivity of aerodynamic power to the blade pitch angle, which can be determined by wind speed, rotor speed, and blade pitch angle. To simplify the model, it is calculated based on the conditions with rated wind speed, rated rotor speed, and zero-pitch angle, and a dimensionless gain-correction factor GK(q) is introduced by:

GKðqÞ ¼

1 1 þ q=qk

> > > > > > KI ðqÞ ¼ > :

IDr U0 u2fn   GKðqÞ vP NGear  ðq ¼ 0Þ vq

(21)

According to the methodology of the GTC and PIC, the wind turbine control strategy is coupled into the LES model. 3. Study site and simulation details 3.1. Study site and data

ðt

Dq ¼ KP NGear DU þ KI NGear DUdt

8 2IDr U0 zf ufn   GKðqÞ > K ðqÞ ¼ > > P vP > > NGear  ðq ¼ 0Þ > > < vq

(20)

where qk is the corrected blade pitch angle at which the pitch angle sensitivity is two times of the pitch angle in the rated operating condition. Ultimately, the exact gains of KP and KI varied with the pitch angle can be calculated by:

To make the model validation and application, a realistic wind farm (Muren) located in the northeastern Inner Mongolia of China, including 33 wind turbines, is selected and simulated using the MNF. The topographical features and the wind turbine layout are displayed in Fig. 3. The altitude in the north is higher than that in the south, and its difference is about 100 m. The WFP model is updated by the technical parameters of a GW-77/1500 wind turbine, as listed in Table 2. The wind turbines are labeled as WT-1, WT-2, …, and WT-33. Besides, the observations including the 15min averaged nacelle wind speed at each wind turbine and the power output during 3 days (3rd to 5th in July 2016) are available. Meanwhile, the observed data have passed the quality control with effective value test, extreme value test, and high consistency test. 3.2. Model configuration Firstly, the mesoscale simulations are initialized by the NCEPFNL data with 6-h intervals. The geographical data resolution is modis_30s þ 30s. A four-dimensional data assimilation system and the grid nudging scheme are used. To obtain the preprocessing results with the highest resolution possible, the model configuration validated by our previous study [47] is applied in this mesoscale simulation as well. For instance, three nested domain is used to the horizontal grid and two-way nesting between the domains is adopted, as shown in Fig. 4, the horizontal resolution is 200 m for the d03 domain, and 81 levels are set in the vertical direction. To diminish the sensitivity of the scheme to the time step and increase the accuracy of the vertical profile, the WSM-6 is applied to the simple microphysical process [48]. The Goddard for short wave radiation scheme and RRTM for long wave radiation schemes are adopted in the radiation process [49,50]. The Grell-Freitas

Q. Wang et al. / Energy 203 (2020) 117913

7

Fig. 3. The terrain of the design domain and the wind farm layout for the microscale simulation.

Table 2 Technical parameters of the wind turbine. Parameter

Value

Parameter

Value

Rated power Rotor diameter Hub height

1.5 MW 77 m 70 m

Rated speed Cut-in/out speed Rated rotor speed

11 m$s1 3/25 m$s1 17.3 r$min1

cumulus parameterization scheme is applied in the thermal diffusion scheme [51]. In addition, the sensitivity of the WRF model to these schemes has been verified and analyzed by Wang et al. [52], showing the atmospheric physical processes can be captured well. Xia et al. [53] proved that MYNN-2.5 [54] is the only PBL scheme that can accurately resolve the turbulence energy caused by wind farms. Thus, the MYNN-2.5 scheme is selected, which can ensure the effective implementation of the WFP [15]. The physical scheme options in the WRF model are summarized in Table 3. In addition, the period from 3rd to 5th in July 2016 is simulated, and the first 24-

h is treated as spin-time. The hourly-average output data are saved. After getting the coarse wind field data in mesoscale considering the interaction between the wind farm and ABL, a fine wind profile with turbulence needs to be established to perform the CFD simulations using the LES model. Typically, a neutral atmosphere is taken into account. Thus, the nocturnal hourly average wind data during 2200LT ~ 2300LT (LT ¼ UTC þ 8 h) on July 4, 2016 from the preprocessor mesoscale simulation is selected and evaluated. Different from the traditional approaches using the conventional logarithmic and power wind profile [55], the cubic polynomial curve is better fitting from the mesoscale wind data, which could serve as a mean inflow wind profile for microscale simulations, as displayed in Fig. 5. Moreover, a turbulence fluctuation is further imposed on this mean inflow wind profile by using the KS method [31], as described in section 2.3. The key parameters of the KS method can be determined based on the wind resource near the Muren wind farm, such as the mean wind speed at 10 m v10,ave is 2.38 m/s and the surface roughness coefficient K is 0.0063 [38].

Fig. 4. The terrain of the three nested domains for the preprocessor mesoscale simulation.

Table 3 Physical scheme options in the mesoscale model. Physical scheme

Option

Physical scheme

Option

Microphysics Short/Long-wave radiation Planetary boundary layer

WSM-6 Goddard/RRTM MYNN-2.5

Cumulus Land-surface model Wind farm parameterization

Grell-Devenyi ensemble United Noah scheme Fitch’s model

8

Q. Wang et al. / Energy 203 (2020) 117913

4. Results and discussion 4.1. Model validation

Fig. 5. The wind profile for mean inflow wind speed.

Thus, a real-like ABL condition is prepared and it serves as the initial inflow condition by compiling the user-defined functions (UDF) in the LES model. Finally, the LES simulations are accomplished based on the open-source platform. The geometric computational domain of the Muren wind farm is established, as displayed in Fig. 3. The domain is 10 km (north-south)  7 km (east-west)  0.4 km (bottom-top). The grid size near the wind turbine blade is 2 m, following the study of Fleming et al. [56], leading to a total grid of 18.54 million. Besides, the turbine induced body force is modeled using the ALM method [25]. The basic parameters in ALM can be determined on the basis of the aerodynamic parameters of the wind turbine GW77/1500 and the criterion given by Troldborg et al. [44]. For example, the actuator point number of the blade, nacelle, and tower is 40, 10, and 80, respectively; the projection width is twice the grid size near the actuator point, and its body force at each point projects uniformly over a disk of some radius but vanish off as a Gaussian function in the directions normal and tangential to the disk. In this case, the relative piecewise generator speed in each stage of these five regions in the pitch and torque controllers are determined by the wind turbine based on the criterion reported by Quallen et al. [33], as labeled in Fig. 2. As mentioned above, this LES simulation is initialized by the prepared cubic polynomial inflow boundary condition. In addition, to meet the CFL law [24] and ensure a fully-developed flow, the time step is 0.025 s and the total computational time is 2700 s. The configuration of the microscale model is summarized in Table 4.

Based on the well-defined configurations of the MNF, the model validation is carried out with the realistic onshore wind farm. Here, both flow property and performance of the wind farm in complex terrain are validated sufficiently. However, it is too time-consuming for the LES simulation to reach a time scale comparable to that of the observation scale. Thus, the 15-min averaged (last 900 s in the total iteration length of 2700 s) data from simulation and observation are compared. In this study, the nacelle velocity and the generator power of each wind turbine simulated by the MNF are compared with the measured data provided by the wind power manufacturer. Fig. 6a shows the scatters of nacelle velocity for each turbine between the results from the simulation and the observation. The coefficient (R2) is 0.9348 and the relative error ε is 3.1%. It can be seen that the simulated nacelle velocity of most wind turbines is slightly overestimated. This is mainly attributed to factors like the high turbulence generated near the nacelle, the error from the ALM method, and the system errors of LES modeling [16]. Comparing with the previous studies by Makridis et al. [57] and Yan et al. [58], a higher linear correlation and lower relative error have been obtained in the present work, confirming that the MNF is better to capture the flow feature of the whole wind farm in complex terrain. In addition, Fig. 6b displays the 15-min averaged velocity contour at the hub height. It shows the wake meandering for each wind turbine and the serious interferences for several wind turbines, such as the speed deficit at downstream of WT-24 and the wake superposition effect at WT-6 and WT-7. Meanwhile, the wind speed enhancement and blockage caused by the wavy terrain are also notably depicted. Further validation is made against the generator power of each wind turbine. Here, to validate the coupled wind turbine control strategy in MNF (MNF-CTR hereafter) as well, the simulations using the framework without control strategies (MNF-NCTR) are also conducted. Fig. 7 and Fig. 8 compare the power of each wind turbine simulated by MNF-CTR and MNF-NCTR with the observed data, respectively. It can be seen that the R2 for the current model is 0.9735, which is larger than that of the MNF-NCTR (0.6561). Both the MAE and the RMSE of the current framework are lower, and the relative error ε is 8.31%, showing the accuracy of the MNF-CTR is 6.65% higher than that of the MNF-NCTR. This improvement mainly contributes to the coupled wind turbine control strategy, which could well model the real-time optimal operating process for each wind turbine when suffering the upwind turbine wake effects [33]. In addition, the model properly reproduces the generator power of each wind turbine with an overestimation. This is partly due to the overestimation of wind speed and also attributed to other factors limiting the power generation by the turbines which are not considered here, such as the peak load regulation of power grid in actual wind farm operation [25]. In summary, although the errors between simulated power and observed data can be seen, the discrepancies are acceptable as indicated by a high linear correlation and low relative error. Therefore, both flow property and

Table 4 The main configurations of the microscale model. Item

Content

Item

Content

Domain size/km Grid resolution/m Time step/s Points (Blade/Nacelle/Tower)

7.5  10  0.4 2 0.025 40/10/80

Point-Dist.-Type Force-Projection-Type Control Strategy Inflow boundary condition

Uniform Disk-Gaussian GTC & PIC UDF

Q. Wang et al. / Energy 203 (2020) 117913

9

Fig. 6. (a) The nacelle velocity comparison between the current simulation and the observation and (b) the velocity contour of the wind farm at the hub height.

Fig. 7. Power of the simulated results (the blue is the MNF-NCTR and the red is the current framework) and observed results (black) for each wind turbine in wind farms. The below depicts the relative error bar. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)

into the flow behavior and power generation of the onshore wind farm in complex terrain. 4.2. Flow behavior and performance

Fig. 8. The power output comparison between the simulations and the observation.

performance of the wind farm in complex terrain are wellvalidated, indicating the MNF-CTR is reasonable to gain insight

The performance of the entire wind farm is mainly determined by the flow behavior near each wind turbine [59]. Based on the averaged flow characteristics of the entire wind farm (Fig. 6b), it can be seen that three adjacent wind turbines (WT-5, WT-6, and WT-7) in line with a typical axial distance of 5D [60] undergo the strongest wake effects. This is a typical local layout scenario in the operation of wind farms in reality, which should be examined. Fig. 9 shows the wake interference property of these wind turbines in three views. Generally, the downwind WT-6 and WT-7 undergo upwind turbine wake effects (Fig. 9a), and the wakes superimpose downwind of WT-6 and meander to WT-7. Fig. 9b presents a significant asymmetry wind speed loss in the vertical direction, i.e., the closer to the ground, the stronger the wind deficit. It can be found in Fig. 9c as well. This mainly contributes to wind shear and terrain effects. In addition, as the evolution of the wake, the asymmetry of the wake becomes more significant in the downstream of WT-6. It can be seen that the #3 blade undergoes the strongest wake, which would induce the non-equilibrium load on the wind turbine, resulting in the lower power output [61]. Furthermore, the performance of these wind turbines is

10

Q. Wang et al. / Energy 203 (2020) 117913

Fig. 9. The wake property of the wind turbines WT-5, WT-6, and WT-7: (a) top view at the hub height, (b) the vertical slice along the stream-wise, and (c) the vertical slices at spanwise.

Table 5 Operation parameters of the selected wind turbines. Wind turbine

U/r$min1

Ta/kN

Taero/kN$m

P/kW

WT-5 WT-6 WT-7

15.1 12.6 12.4

147.3 85.8 80.8

631.9 334.2 314.2

1116.0 491.0 455.2

examined, and the output parameters like rotor speed U, axial thrust Ta, rotor torque Taero, and power output P are listed in Table 5. The sharp reduction of operational performance occurs on the downstream wind turbines. Taking the power as an example, WT-6 and WT-7 suffer the power loss of 56.0% and 59.2%, respectively. This is largely determined by the inflow speed. Quantitatively, the inflow profile of these wind turbines is displayed in Fig. 10. Here, the above ground level is normalized by the rotor diameter (z/D), and the shadow stands for the rotor area. It can be seen that the maximum wind deficit reaches 4 m/s near the hub-heights of WT-6 and WT-7. Accordingly, this typical local layout scenario is unfavorable to the output of the whole wind farm. A good suggestion is to control the yaw of the upstream turbine (WT-5) to make its wake deflect, reducing the impact on the downstream [62]. In addition, a speed-up zone in the range from land surface to the lowest edge of the rotor for these wind turbines, especially for the vicinity of WT-7. It means this acceleration effect caused by terrain lift may compensate for the wake effects caused by WT-6. More evidence about the mechanism of the influence of topography on the operation process of the wind turbine needs to be detected.

4.3. Operating property in complex terrain Recently, wind farm development in complex terrain becomes more and more extensive [63]. As mentioned above, there exists a speed-up zone over the wavy terrain and its impact on the operation process of the wind turbine is essential to be explored. Moreover, because the magnitude of the wind in nature always

Fig. 10. Wind profiles at the position of the three wind turbines.

changes, it is significant to gain insight into the operating characteristics of the wind turbine over wavy terrain under a dynamic inflow condition. To assess the impact of terrain on the operating process of wind turbines in complex terrain, a dynamic inflow wind with a sinusoidal-temporal variation is established based on the mentioned cubic polynomial wind profile. This profile can be expressed as a function of height z and time t:

Q. Wang et al. / Energy 203 (2020) 117913

Uin ðz; tÞ ¼ Upoly ðzÞ,ðcv þ jsinðpt = TÞjÞ

(22)

where Upoly is the cubic polynomial fitting wind speed, z is the height above ground level, cv ¼ 0.3 is the inflow correction constant, and T is an inflow fluctuation period (T ¼ 300 s). Other configurations are the same as mentioned above, and the total iteration length is four periods (4T). Here, three wind turbines located in three typical locations of hills were selected to analyze the operation characteristics and influence mechanism of wind turbines in different locations of complex terrain. The selected wind turbines are WT-8, WT-12, and WT-10, respectively defined as the upwind wind turbine (UWT), top wind turbine (TWT), and downwind wind turbine (DWT), as illustrated in the left panel of Fig. 11. Meanwhile, the averaged (last period in 4T) wind contour is displayed, where the Uhub represents the axial velocity at hub-height and CUrotorD represents the spatial average velocity for the rotor area. It shows that the inflow velocity for TWT is the largest, and the wake recovery of TWT is the fastest. This is followed by UWT, and DWT is the slowest. This discrepancy will further change the operating process and performance, such as the variations of the rotor speed, pitch angle, rotor/generator torques, driving force, and generator power. Taking the last period as an example, the variations of rotor speed and pitch angle with the change of inflow wind speed for UWT, TWT, and DWT are depicted in Fig. 12a. The time when the inflow wind speed is greater than the rated value is recorded as trate-in and lower than the rated value is recorded as trate-out. In theory, the PIC strategy begins to work at trate-in and ends at trate-out. However, it is noted that the initial time and terminal time of PIC for the wind turbines over complex terrain are advanced or delayed due to the terrain effect. For TWT, as the inflow rises, the rotor speed increases sharply and then reaches the rated value, which is apparently earlier than UWT and DWT. Then, the rotor speed is limited and slows down slightly until the inflow speed is lower than the rated value. After the trate-out, the rotor speed is reduced sharply and the PIC stops working (the pitch angle plummeted to zero). Meanwhile, the rotor speed of TWT diminishes more slowly than UWT and DWT. It concludes that the performance of the UWT or

11

DWT is weakened due to the turbulence caused by the hill. In contrast, the TWT draws support from the acceleration effect, prolongs the optimal operation time, and promotes the operating performance. It is acknowledged that the variability of the rotor speed is regulated by the GTC strategy. The variation of aerodynamic and generator torques is illustrated in Fig. 12b. Taking the TWT as an example, at the first stage, the torques increase sharply as the inflow increases and reaches the maximum value before the trate-in. The magnitude of the rotor aerodynamic torque is greater than the generator torque, according to Eq. (17), thus the rotor speed increases. When the wind speed is larger than the rated value, both two torques quickly drop to be nearly equal and fluctuate within a slight range until it is lower than the rated value again. Subsequently, the torques decrease. While the generator torque is larger than the aerodynamic torque, resulting in the reduction of the rotor speed. Similar trends in the GTC process can be found for UWT and DWT. Consequently, the parameter changes as the inflow varies, such as the axial force Ta (Fig. 13a) and the generator power P (Fig. 13b). Notably, more power generation is produced by TWT than that of UWT and DWT. A fundamental reason is that the speed-up zone occurs at the top of the hill, the blocking zone generates at the upwind, and the vortex cavity with high turbulence induces at the downwind [64]. Further, the operating difference between these wind turbines is quantified. The change of time length proportion when the generator power cut in or cut out the rated-state (DTPs) is defined as:

DTPs ¼

trates  ts T

(23)

where ts represents the time when the generator power cut in or cut out the rated value, and trate-s can be replaced by trate-in or trateout. The positive value of DTPin or DTPout means reaching the ratedstate in advance, contrarily, the negative indicates delay. Thus, DTPtot ¼ DTPin  DTPout stands for the change of the rated-state operating time in one period, as shown in Fig. 14. It shows that the time when the generator power cut in the

Fig. 11. The flow behavior for the wind turbines on different positions of the hill (a) UWT, (b) TWT, and (c) DWT.

12

Q. Wang et al. / Energy 203 (2020) 117913

Fig. 12. The variations of the (a) rotor speed and pitch angle and (b) rotor and generator torque for the wind turbines of UWT, TWT, and DWT.

Fig. 13. The variations of the performance parameters such as (a) rotor axial force and (b) generator power for three wind turbines of UWT, TWT, and DWT.

the total operating time length proportion under the rated-state for TWT increases by 5.5%. It means that the TWT can cut in the ratedstate earlier and cut out later, so that the working hours under the rated condition in the unit cycle is increased, thereby improving the output. But that for UWT and DWT decreases by 1.0% and 3.0%, respectively. Overall, the DTPtot for these three wind turbines prolongs by 1.0%. It indicates that the generation of the entire wind farm could be optimized through making full use of terrain acceleration effects in complex terrain. 5. Conclusions

Fig. 14. The change of the time length proportion when the generator power cut-in or cut-out the rated-state (DTPin or DTPout) and the total change of the time length proportion (DTPtot) for UWT, TWT, and DWT.

rated-state advances 3.1% and the time length when the generator power cut out the rated-state delays 2.4% for TWT. In other words,

A high-fidelity multiscale numerical framework to assess the performance of the wind farm in complex terrain in the real atmosphere, combining the WRF model coupled with a wind farm parameterization and a LES model embedded with the wind turbine control strategy is developed. The model validation and application are carried out at a realistic onshore wind farm in complex terrain located in Inner Mongolia of China. The nacelle velocity and the generator power of each wind turbine simulated by the MNF are compared well with the measured data, e.g., the coefficient R2 and the relative error ε of nacelle velocity are 0.9348 and 3.1%, respectively. A typical local layout scenario in the operation of wind farms in reality, e.g., three adjacent wind turbines in line with a typical axial

Q. Wang et al. / Energy 203 (2020) 117913

distance of 5D, undergoes a strong wake effect. The downwind wind turbines suffer a power loss of 56.0% and 59.2%, respectively. Accordingly, this typical local layout scenario is unfavorable to the output of the whole wind farm. A good suggestion is to control the yaw of the upstream turbine to make its wake deflect, mitigating the impact of non-equilibrium load on the downstream wind turbines. In addition, the impact of terrain on the operating property under dynamic inflow condition is quantified. The wind turbine at the top of the hill can cut in the rated-state earlier and cut out later, increasing the working time under the rated state in the unit cycle, thereby improving the power output of 5.5%. Overall, the total rated-state operating time length for these wind turbines is extended by 1.0%, showing it can improve the power generation of the entire wind farm via making full use of the terrain acceleration effects in complex terrain. This multiscale numerical framework has been well-validated and it can be used to assess the flow behavior and performance of a realistic onshore wind farm in complex terrain, which facilitates the siting of wind farms, especially for the distributed wind power. It can also serve as a tool for evaluating the operating process under dynamic inflow conditions. However, how to realize the online coupling between the mesoscale and the microscale models remains to be a challenge in this field, and needs more efforts in future research. CRediT authorship contribution statement Qiang Wang: Writing - original draft, Visualization, Investigation. Kun Luo: Conceptualization, Methodology, Writing - review & editing. Renyu Yuan: Software, Validation. Shuai Wang: Writing review & editing. Jianren Fan: Supervision. Kefa Cen: Supervision. Acknowledgments This work was supported by the National Youth Top-notch Talent Support Program and the National Natural Science Foundation of China (Grants No. 51766014). We are grateful to that and the members who are not listed as co-authors for their helpful discussions and comments. References €hme GS, Fadigas EA, Gimenes ALV, et al. Wake effect measurement in [1] Bo complex terrain - a case study in Brazilian wind farms. Energy 2018;161: 277e83. [2] Ren G, Liu J, Wan J, et al. Overview of wind power intermittency: impacts, measurements, and mitigation solutions. Appl Energy 2017;204:47e65. [3] Willis D, Niezrecki C, Kuchma D, et al. Wind energy research: state-of-the-art and future research directions. Renew Energy 2018;125:133e54. [4] Blocken B, Hout AVD, Dekker J, et al. CFD simulation of wind flow over natural complex terrain: case study with validation by field measurements for Ria de Ferrol, Galicia, Spain. J Wind Eng Ind Aerod 2015;147:43e57. [5] Risan A, Lund J, Chang C-Y, et al. Wind in complex terrain-lidar measurements for evaluation of CFD simulations. Rem Sens 2018;10(59):1e18. [6] Subramanian B, Chokani N, Abhari RS. Aerodynamics of wind turbine wakes in flat and complex terrains. Renew Energy 2016;85:454e63. [7] Bai C-J, Wang W-C. Review of computational and experimental approaches to analysis of aerodynamic performance in horizontal-axis wind turbines (HAWTs). Renew Sustain Energy Rev 2016;63:506e19. [8] Siuta D, West G, Stull R. WRF hub-height wind forecast sensitivity to PBL scheme, grid length, and initial condition choice in complex terrain. Weather Forecast 2017;32(2):493e509. [9] Mughal MO, Lynch M, Yu F, et al. Wind modelling, validation and sensitivity study using Weather Research and Forecasting model in complex terrain. Environ Model Software 2017;90:107e25. [10] Dzebre DEK, Adaramola MS. A preliminary sensitivity study of Planetary Boundary Layer parameterization schemes in the weather research and forecasting model to surface winds in coastal Ghana. Renew Energy 2020;146: 66e86. mez-Gesteira M, et al. A sensitivity study of the WRF [11] Carvalho D, Rocha A, Go model in wind simulation for an area of high wind energy. Environ Model Software 2012;33:23e34.

13

 mez-Gesteira M, et al. Sensitivity of the WRF model [12] Carvalho D, Rocha A, Go wind simulation and wind energy production estimates to planetary boundary layer parameterizations for onshore and offshore areas in the Iberian Peninsula. Appl Energy 2014;135:234e46. [13] Baidya Roy S. Can large wind farms affect local meteorology? J Geophys Res 2004;109(D19101):1e6. [14] Blahak UGB, Meis J. A simple parameterization of drag forces induced by large wind farms for numerical weather prediction models. 2010. [15] Fitch AC, Olson JB, Lundquist JK. Parameterization of wind farms in climate models. J Clim 2013;26(17):6439e58. -Agel F. A new wind-farm parameterization for large-scale [16] Abkar M, Porte atmospheric models. J Renew Sustain Energy 2015;7(13121):1e12. nez PA, Navarro J, Palomares AM, et al. Mesoscale modeling of offshore [17] Jime wind turbine wakes at the wind farm resolving scale: a composite-based analysis with the Weather Research and Forecasting model over Horns. Rev. Wind Energy. 2015;18(3):559e66. [18] Lee JCY, Lundquist JK. Observing and simulating sind-turbine wakes during the evening transition. Bound-lay Meteorol. 2017;164(3):449e74. [19] Lee JCY, Lundquist JK. Evaluation of the wind farm parameterization in the Weather Research and Forecasting model (version 3.8.1) with meteorological and turbine power data. Geosci Model Dev (GMD) 2017;10(11):4229e44. [20] Mangara RJ, Guo Z, Li S. Performance of the wind farm parameterization scheme coupled with the Weather Research and Forecasting model under multiple resolution regimes for simulating an onshore wind farm. Adv Atmos Sci 2018;36(2):119e32. [21] Toja-Silva F, Kono T, Peralta C, et al. A review of computational fluid dynamics (CFD) simulations of the wind flow around buildings for urban wind energy exploitation. J Wind Eng Ind Aerod 2018;180:66e87. [22] Mehta D, van Zuijlen AH, Koren B, et al. Large Eddy Simulation of wind farm aerodynamics: a review. J Wind Eng Ind Aerod 2014;133:1e17. -Agel F, Bastankhah M, Shamsoddin S. Wind-turbine and wind-farm [23] Porte flows: a review. Bound-lay Meteorol 2019:1e59. [24] Martinez L, Leonardi S, Churchfield M, et al. A comparison of actuator disk and actuator line wind turbine models and best practices for their use. In: 50th AIAA aerospace sciences meeting including the new horizons forum and aerospace exposition; 2012. [25] Stevens RJAM, Martínez-Tossas LA, Meneveau C. Comparison of wind farm large eddy simulations using actuator disk and actuator line models with wind tunnel experiments. Renew Energy 2018;116:470e8. [26] Churchfield MJ, Sang L, Moriarty PJ. Adding complex terrain and stable atmospheric condition capability to the OpenFOAM-based flow solver of the simulator for on/offshore wind farm applications (SOWFA). Golden, CO (United States): National Renewable Energy Lab; 2013. [27] Sharma V, Cortina G, Margairaz F, et al. Evolution of flow characteristics through finite-sized wind farms and influence of turbine arrangement. Renew Energy 2018;115:1196e208.  M, Mouslim H, et al. From meso-scale to micro scale LES [28] Medjroubi W, Mache modelling: application by a wake effect study for an offshore wind farm. ITM Web of Conferences 2014;2:01004. [29] Carvalho D, Rocha A, Santos CS, et al. Wind resource modelling in complex terrain using different mesoscaleemicroscale coupling techniques. Appl Energy 2013;108:493e504. [30] Gopalan H, Gundling C, Brown K, et al. A coupled mesoscaleemicroscale framework for wind resource estimation and farm aerodynamics. J Wind Eng Ind Aerod 2014;132:13e26. [31] Meyer DW, Eggersdorfer ML. Simulating particle collisions in homogeneous turbulence with kinematic simulation e a validation study. Colloids & Surfaces A Physicochemical & Engineering Aspects 2014;454(454):57e64. [32] Fung JCH, Hunt JCR, Malik NA, et al. Kinematic simulation of homogeneous turbulence by unsteady random Fourier modes. J Fluid Mech 2006;236: 281e318. [33] Quallen S, Xing T. CFD simulation of a floating offshore wind turbine system using a variable-speed generator-torque controller. Renew Energy 2016;97: 230e42. [34] Skamarock WC, Klemp JB, Dudhia J, et al. A description of the Advanced Research WRF version 3.0 (NCAR TECHNICAL NOTE). National Center for Atmospheric Research; 2008. [35] Fitch AC, Olson JB, Lundquist JK, et al. Local and mesoscale impacts of wind farms as parameterized in a mesoscale NWP model. Mon Weather Rev 2012;140(9):3017e38. [36] Meyer DW, Eggersdorfer ML. Simulating particle collisions in homogeneous turbulence with kinematic simulation-A validation study. Colloid Surface Physicochem Eng Aspect 2014;454:57e64. [37] Ducasse L, Pumir AJPRE. Inertial particle collisions in turbulent synthetic flows: quantifying the sling effect 2009;80(6):066312. [38] Pan T. Numerical simulation of atmospheric boundary layer wind field by OpenFOAM. Chongqing: Chongqing University; 2015. -Agel F. Modeling turbine wakes and power losses within a [39] Wu Y-T, Porte wind farm using LES: an application to the Horns Rev offshore wind farm. Renew Energy 2015;75:945e55. [40] Tominaga Y, Mochida A, Murakami S, et al. Comparison of various revised keε models and LES applied to flow around a high-rise building model with 1:1:2 shape placed within the surface boundary layer. J Wind Eng Ind Aerod 2008;96(4):389e411. [41] Tominaga Y, Stathopoulos T. Numerical simulation of dispersion around an

14

[42]

[43] [44]

[45]

[46]

[47]

[48] [49]

[50] [51] [52]

Q. Wang et al. / Energy 203 (2020) 117913 isolated cubic building: model evaluation of RANS and LES. Build Environ 2010;45(10):2231e9. Tominaga Y, Stathopoulos T. CFD simulation of near-field pollutant dispersion in the urban environment: a review of current modeling techniques. Atmos Environ 2013;79:716e30. Sorensen JN, Shen WZ. Numerical modeling of wind turbine wakes. J Fluid Eng 2002;124(2):393. Troldborg N, Larsen GC, Madsen HA, et al. Numerical simulations of wake interaction between two wind turbines at various inflow conditions. Wind Energy 2011;14(7):859e76. Li Qa, Kamada Y, Maeda T, et al. Fundamental study on aerodynamic force of floating offshore wind turbine with cyclic pitch mechanism. Energy 2016;99: 20e31. J. Jonkman Sb, W. Musial, and G. Scott. Definition of a 5-MW reference wind turbine for offshore system development. Prepared under Task No WER533012009. Yuan R, Ji W, Luo K, et al. Coupled wind farm parameterization with a mesoscale model for simulations of an onshore wind farm. Appl Energy 2017;206:113e25. Hong S-Y, Noh Y, JJMwr Dudhia. A new vertical diffusion package with an explicit treatment of entrainment processes 2006;134(9):2318e41. Mlawer EJ, Taubman SJ, Brown PD, et al. Radiative transfer for inhomogeneous atmospheres: RRTM, a validated correlated-k model for the longwave 1997;102(D14):16663e82. Fouquart Y, Bonnel B. Ramaswamy VJJoGRA. Intercomparing shortwave radiation codes for climate studies 1991;96(D5):8955e68. ve nyi DJGRL. A generalized approach to parameterizing convecGrell GA, De tion combining ensemble and data assimilation techniques 2002;29(14). 38-4. Wang C, Jin S. Error features and their possible causes in simulated low-level winds by WRF at a wind farm. Wind Energy; 2013. p. 1315e25.

[53] Xia G, Cervarich MC, Roy SB, et al. Simulating impacts of real-world wind farms on land surface temperature using the WRF model: validation with observations. Mon Weather Rev 2017;145(12):4813e36. [54] Nakanishi M, HJJotMSoJSI Niino. Development of an improved turbulence closure model for the atmospheric boundary layer 2009;87(5):895e912. [55] Shu ZR, Li QS, He YC, et al. Observational study of veering wind by Doppler wind profiler and surface weather station. J Wind Eng Ind Aerod 2018;178: 18e25. [56] Fleming PA, Gebraad PMO, Lee S, et al. Evaluating techniques for redirecting turbine wakes using SOWFA. Renew Energy 2014;70:211e8. [57] Makridis A, Chick J. Validation of a CFD model of wind turbine wakes with terrain effects. J Wind Eng Ind Aerod 2013;123:12e29. [58] Yan BW, Li QS. Coupled on-site measurement/CFD based approach for highresolution wind resource assessment over complex terrains. Energy Convers Manag 2016;117:351e66. [59] Meneveau C. Big wind power: seven questions for turbulence research. J Turbul 2019;20(1):2e20. [60] Kuo JYJ, Romero DA, Beck JC, et al. Wind farm layout optimization on complex terrains e integrating a CFD wake model with mixed-integer programming. Appl Energy 2016;178:404e14. [61] Ciri U, Rotea MA, Leonardi S. Effect of the turbine scale on yaw control. Wind Energy 2018;21(12):1395e405. [62] Gebraad PMO, Teeuwisse FW, van Wingerden JW, et al. Wind plant power optimization through yaw control using a parametric model for wake effectsa CFD simulation study. Wind Energy 2016;19(1):95e114. [63] Yang X, Pakula M, Sotiropoulos F. Large-eddy simulation of a utility-scale wind farm in complex terrain. Appl Energy 2018;229:767e77. [64] Yan S, Shi S, Chen X, et al. Numerical simulations of flow interactions between steep hill terrain and large scale wind turbine. Energy 2018;151:740e7.