Analyzing complex wake-terrain interactions and its implications on wind-farm performance.

Rotating wind turbine blades generate complex wakes involving vortices (helical tip-vortex, root-vortex etc.).These wakes are regions of high velocity deficits and high turbulence intensities and they tend to degrade the performance of down-stream turbines. Hence, a conservative inter-turbine distance of up-to 10 times turbine diameter (10D) is sometimes used in wind-farm layout (particularly in cases of flat terrain). This ensures that wake-effects will not reduce the overall wind-farm performance, but this leads to larger land footprint for establishing a wind-farm. In-case of complex-terrain, within a short distance (say 10D) itself, the nearby terrain can rise in altitude and be high enough to influence the wake dynamics. This wake-terrain interaction can happen either (a) indirectly, through an interaction of wake (both near tip vortex and far wake large-scale vortex) with terrain induced turbulence (especially, smaller eddies generated by small ridges within the terrain) or (b) directly, by obstructing the wake-region partially or fully in its flow-path. Hence, enhanced understanding of wake- development due to wake-terrain interaction will help in wind farm design. To this end the current study involves: (1) understanding the numerics for successful simulation of vortices, (2) understanding fundamental vortex-terrain interaction mechanism through studies devoted to interaction of a single vortex with different terrains, (3) relating influence of vortex-terrain interactions to performance of a wind-farm by studying a multi-turbine wind-farm layout under different terrains. The results on interaction of terrain and vortex has shown a much faster decay of vortex for complex terrain compared to a flatter-terrain. The potential reasons identified explaining the observation are (a) formation of secondary vortices in flow and its interaction with the primary vortex and (b) enhanced vorticity diffusion due to increased terrain-induced turbulence. The implications of this vortex-terrain interactions on wind-farm performance is observed by comparing two LES simulations of a multi-turbine wind-farm layout (in real actual complex terrain and a made-up flat terrain scenario) with the observed annual power data at the actual wind-farm. The comparison reveals drop in power production due to terrain and wake effects for flatter terrain case. The insights from this study can serve as a step towards enhancing wake-dissipation through either artificial obstruction or artificial terrain modifications.


Introduction and work objective 1.Introduction
Wind-flows in highly complex terrain leads to terrain induced turbulence. The spatio-temporal scales of flow-structures generated in terrain-induced flows are dependent upon terrain features (their form, slope, height and orientation) and incoming flow conditions (wind speed, direction and stratification). These flow-structures influence wake development behind a wind-turbine located in complex terrain. Terrain features that can influence wind-flow can be identified through the study of altitude contour level, and comprises of major terrain features (like, hill, ridge, valley, saddle, depression) and minor terrain features (like spur and cliff). Sharp edges in these features can make the wind-flow to separate and form eddies. The air-wake can establish on the leeway side and lead to airflow reversal. The interaction between wind-turbine induced wakes and terrain induced turbulence has not yet been studied in detail. It is expected that eddy generated due to flow-separation across terrain can interact with the wake region leading to faster wake recovery. Wind-turbine wake region has some distinct flow features -categorized as near wake region (region influenced by blade design that creates distinct tip-vortex near blade-tip, which separates inner wake region from outer free flow) and far-wake region (where large scale turbulence exist owing to wake-expansion and interaction/break-down of helical tipvortices). The near-wake helical tip vortices over a period of many rotor diameters interact and breakdown leading to a far-wake region. Such tip-vortices (and other vortex structures) have been represented by simple wake vortex models (like, Rankine model) [1,2,3] and same has been done in this work. The terrain can influence any part of wake region either by directly obstructing it fully or partially and indirectly by terrain induced turbulence. Keeping in mind the scarcity of research in this area, the objective of the work is to:

Objectives
(i) To get a better insight into wake-terrain interaction and the associated decay of vortex.
(ii) Understanding influence of numerics (meshing, discretization and temporal schemes) in capturing wake-dynamics (near and far-wake) and avoiding numerical diffusion. (iii) Understanding influence of wake-terrain interactions on a multi-turbine wind-farm. (iv) Serve as an initial step towards (a)devising alternative ways for initiating early wakedissipation using obstacles and (b)establishing complex terrain-specific guideline for turbine layout.

Methodology : validation and simulation details
The Large Eddy Simulation (LES) turbulence model used in our earlier work [4] has been chosen to understand the behavior of wake-vortex-terrain interactions. In this work, we isolate the sole influence of terrains in influencing wake dynamics (especially near-wake tip vortex). The nearwake tip vortex has been represented by simple wake-model (like, Rankine as done in some notable work like [1,2,3]). The approach involves studies of interaction of such vortices in wake (represented by an Rankine vortex model) with terrain of different complexities (like, flat-terrain and idealized hill) as shown in figure 1. The challenging part of simulation is to avoid numerical diffusion of vorticity. Hence, proper validation is required to check that numerical schemes and mesh-size do not enhance diffusion. This is done by validating for an inviscid Lamb-Oseen wakevortex behavior (as described in section 2.1). The information from this validation study has been incorporated in conducting LES simulation of wake-vortex-terrain interaction (as in section 2.2). The OpenFOAM-solver is known to be validated for simulation over bumps in a channel (which is similar to Hunt-Hill geometry case here), so we focus on vortex validation only. The equations describing the LES model (a one-equation sub-grid scale turbulent kinetic energy LES model) are given next followed by a description of the validation work. The equations for LES are derived by applying filtering operator to the Navier-Stokes equations. The equation set is represented by Equations 1-2, where velocity u represents the filtered (or resolved) velocity. The three equations below represent a filtered mass continuity equation (Equation 1) and filtered momentum transport equation (Equation 2).
where ρ is the density. where, operator D Dt refers to total derivative, operator ∇ refers to computing gradient, operator ∇· refers to computing divergence, p is pressure, t is time, f refers to external forces arising from actuator line model, R is referred to turbulent stresses and arises owing to averaging procedure. Components of R can be computed as where subscripts i, j refers to components of vector, k is turbulent kinetic energy and ν T is turbulent diffusivity. The computation of eddy viscosity is needed for the closure of equation set and it is obtained from equation 3.
where C k is a constant, ∆ is filter width and is obtained explicitly as cube root of volume of mesh cell. The k sgs is sub-grid scale turbulent kinetic energy and is obtained by solving its transport equation (4). The value of C k used in this study is 0.094.
where,  Figure 1. Case set-up at initial time t = 0 for studying interaction of wake-vortex with flat terrain, hill geometry and scaled hill geometry.

Validation and verification through inviscid Lamb-Oseen wake-vortex behavior
Validation involves numerical study of dynamic wake-vortex behavior in inviscid conditions (zero-viscosity case). Under inviscid condition, the wake-vortex will not decay in time and it will preserve its shape and structure (as seen by vorticity contours). The equation for Lamb-Oseen vortex and its analytical solution is given below. In terms of the vorticity ω, and flow velocity in the circumferential θ-direction V θ , the initial condition is given by and the analytical solution of the Navier-Stokes equations is given by Here, we use initial circulation (Γ 0 = 1), r the radial distance from vortex center and core radius (r c = d c /2 = 0.15). In the absence of viscosity the Lamb-Oseen vortex presents a steady state solution. This validation study helps to understand influence of mesh-size and numerical scheme in avoiding numerical diffusion of wake-vortex. The results are validated with an analytical solution and is verified with a vortex-particle method (VPM)( [5] for time up-to 1 s (as shown in figure 2). The 2D laminar simulation involves vortex size of diameter d c = 0.3m in computational domain of range [−1, 1]m in X-Y direction. Grid-independent solution is achieved when about 130 numerical cells are used across the core radius diameter (d c /δx = 130). The vorticity contour and error contour result in figure 2 indicate that the grid-size and numerical schemes (second order linear convection scheme) are sufficient to avoid high levels of numerical diffusion and provide accurate solutions. The CFD solutions are comparable with analytical model and VPM results. The result from grid-independence reveals that when we use a gridsize one-third of the grid-independent size, then the error increases marginally. Still, this error for zero-viscosity case (mostly attributed to numerical diffusion) will definitely be much less than the physical diffusion in turbulent flows (where, viscous and turbulent diffusion will be high). Hence, for our wake-vortex-terrain interaction studies involving turbulent flows, we use a grid-size of (d c /δx = 50). Further, figure shows adequacy of mesh for LES Simulations through contour of ratio of sub-grid-scale-kinetic-energy k sgs to total-kinetic-energy k total at initial and different times for the 3 cases. In our LES simulations, the mesh-size determines spatial filtering and establishes cut-off between resolved and unresolved (modeled sub-grid scale) parts of flow. Finer the mesh, more part of flow is resolved and more accurate are the simulation. As per criteria of Pope (2000) [6], for a well-resolved LES, less than 20% of the total kinetic energy should be modeled sub-grid-scale part (i.e. k sgs /k total ratio should be less than 0.2). The total kinetic energy k total comprises of the resolved kinetic energy and the modeled sub-grid-scale turbulent kinetic energy. Figure 3 shows that for flatter-terrain, this criteria is satisfied in all the regions at all times (results shown for time t=0 s and t=10 s) with maximum value of ratio reaching 0.15 at the wall. For the hill cases too, most of the regions have satisfied this criteria including the core wake-vortex region except the region extremely close to terrain (where max value is around 0.35) as seen in Figure 3. The mesh can be considered to be adequately resolved for LES simulation of vortex-terrain interaction.

Wake-Vortex-terrain interaction: CFD set-up and solver details
The approach for studying wake-terrain interaction involves, first, obtaining a pre-cursor LES field for all the three terrain cases. The precursor LES field is then considered as an initial field over which Rankine vortex is super-imposed. The simulations are then conducted to analyze the interaction between vortex and terrain induced flow. The computational domain for LES simulation of terrain-vortex interaction is shown in figure 1. A domain size of 5m × 2m × 2m in stream-wise X-direction, span-wise Y-direction and vertical Z direction respectively is used. Along the streamline a periodic boundary condition is imposed and a pressure gradient of 5 × 10 −5 Kg/ms 2 is used to drive the flow. The bottom surface is modelled as a no slip boundary while for the rest of the surfaces a slip boundary condition is used. The profiles given by Equations 8, 9 and 10 were extruded in the span-wise direction to create the bottom terrain surfaces for the three cases 1,2 and 3. It should be noted that the profile in Case 3 is just a 0.66 times scaled down profile of the one used in case 2.
where H = 0.25m. In Case 2, the profile is similar to the one used by Hunt [7] with ratio of hill height/hill length=0.05. However, the actual hill used in that [7] was obtained by rotating the  Figure 3. Adequacy of mesh for LES Simulations -contour of ratio of sub-grid-scale-kineticenergy to total-kinetic-energy at initial and different times for the 3 cases.
profile given by Equation 9 . The peak of the hill is located at about 1.5m from the upstream periodic boundary (ratio of distance from upstream periodic inlet/hill height = 0.166) and the downstream periodic boundary is located 3.5m away from the hill peak (ratio of distance from downstream periodic outlet /hill height = 0.07). The location of the hill and the domain size is thus sufficient to avoid periodicity of hill generated flow structures. The wake-vortex is superimposed on a precursor LES field for each case at time t = 0. The precursor LES simulation was carried out till the flow appeared to have developed in each case (which took around 10 flow through times). The super-imposed wake-vortex has an overall diameter d c = 0.365m and an inner viscous core of 0.091m. The center of the wake-vortex is located at a distance of about 3.5 times the hill height from the terrain. The mesh size for this domain is about 2 million cells with most cells made of hexahedral elements. The grid-size (δx) is such that d c /δx = 50 (as mentioned is section 2.1) but near the terrain, the grid is made finer. The information about adequacy of mesh is mentioned in section 2. In this work, all the equations (except k and turbulence equations) use second order linear discretization scheme, while the turbulent equations use hybrid linear-upwind convection schemes. A backward temporal discretization scheme is used. Similarly, the diffusion term involving Laplacian operator (the divergence of the gradient) is simplified to computing the gradient of the variable at the face. The gradient term can be split into contribution from the orthogonal part and the non-orthogonal parts, and both these contributions have been accounted for. The influence of terrain on wake-vortex decay is discussed in the results section.
In the rest of the paper the Case 1, 2 and 3 are referred to as flat, Hunt Hill and Scaled Hunt Hill respectively.

Formation of Secondary Vortices and its interaction with the Primary Vortex
The decay of vortex on interaction with terrain is discussed in this section. Figure 4 shows vorticity iso-surfaces depicting vortex-terrain interaction dynamics for the three cases. The isosurfaces at different times (t = 1s, 2s, 3s, 7s and 10s) can be seen in the figure. The results indicate that vortex decay is the fastest for Hunt Hill case (as seen at 10s time) followed by Scaled Hunt Hill and the flat terrain case. The vortex decay in all the 3 cases can be attributed to the formation of secondary vortices (SV) and its interaction with the primary vortex. SV's can be identified in the figure as a ring of vortices enclosing the main vortex and they are seen to be forming as early as time t = 1s for all the cases. The SV's are formed as a result of shear conditions in the flow due to vortex -ambient flow interaction (which are influenced by terrain geometry). The number of SV's are less in case of flat-terrain at t = 1s as compared to the Hunt hill and to the Scaled Hunt hill case. Further, in cases 2 and 3, more SV's are concentrated in region vertically above the hill peak height that the flat-terrain case (as seen at time t = 1s). Further, the magnitude of vorticity in SV's and on the surface of primary vortex is higher in case with Hunt Hill than the scaled Hunt Hill and flat terrain case. The interaction between SV and primary vortex leads to annihilation of vorticity. The mechanism of SV formation and interaction is the reason behind faster decay of primary vortex in case of Scaled Hunt hill case (as seen at 10s time) followed by Scaled Hunt hill and flat terrain. The turbulence generated during the process also contributes to the vorticity diffusion. This can be seen in figure 5, which shows instantaneous total turbulent kinetic energy (resolved + sub-grid-scale) for the 3 cases.  Figure 4. Vorticity iso-surface depicting wake-terrain interaction dynamics for 3 cases -the flat terrain, Hunt hill and Scaled Hunt hill.

Turbulent Stresses
A more detailed analysis can be made from observing the turbulent stress components (as in figure 6) and turbulent stress magnitude (as in figure7). Turbulent stress components reveal a dominant mechanism for vorticity-decay through diffusion in flat-terrain and Hunt Hill case. At t = 1s, a higher value of turbulent kinetic energy (figure 5) and turbulent stresses (figures 6-7) can be observed in regions over the hill location and in the vortex region above the hills as compared to the flat-terrain case. A comparison of turbulent stress components magnitudes for flatter terrain and Hunt hill case (from figure 6) reveals that in both the cases the dominant turbulent diffusion direction is in radial direction (YY, ZZ, YZ), and the magnitudes of turbulent stress components is higher in case of Hunt hunt hill than flatter terrain (thus suggesting a higher turbulent diffusion). These turbulent stresses could also be representing interaction between primary anc secondary vortices apart from playing a role in diffusion of vortices. The magnitude of peak turbulent kinetic energy and peak turbulent stresses goes on reducing with time (1 to 10 s) in all the three cases as the primary vortex decays away. At t = 10s in figure 5, turbulent kinetic energy contour confirms that vortex is almost destroyed in case of full-scale Hunt Hill and it has the fastest decay, while flatter-terrain case has slowest wake decay. Another interesting thing revealed in figure 5 at t = 7 − 10s is a possible interaction between decaying vortex and separating flow from hunt hill case. The overall results indicate the role of terrain / obstruction  Figure 5. Instantaneous total turbulent kinetic energy (resolved + sub-grid-scale) depicting wake-terrain interaction dynamics for 3 cases -the flat terrain, Hunt hill and Scaled Hunt hill.

Flat Terrain
Hunt Hill Figure 6. Resolved turbulent stress components depicting dominant momentum transfer mechanism for wake-decay in flat-terrain and hunt hill case.

Implications for wind-farm performance
Influence of terrain on wake and power production is studied for the real multi-turbine Bessaker wind farm (as shown in figure 9). The observed power data from the actual operating Bessaker wind-farm in complex terrain is used for comparison with the predictions from LES simulations. The computational set-up for LES simulations involves : (a) the industrial Bessaker wind-farm situated in its actual complex terrain and (b) the same Bessaker wind-farm layout in a made-up flat terrain. The details of computaional methodology, i.e. the LES simulation methodology and the Actuator line control methodology (i.e. used to model the turbine behavior) are as described in our previous work [4]. Figure 10 shows the influence of wake and terrain effects through comparison of power production per turbine for the same wind-farm turbine layout in flat terrain and in realistic complex terrain scenarios. The observed power data is perturbed (scaled slightly) so as to maintain the confidentiality of exact power generated (as requested by the industry) but we ensure that the shown power trend vis-a-vis the simulated results (as in figure 10) and the final conclusions are not affected by this scaling.  Figure 9. Bessaker wind farm in complex terrain (left) and the observed annual power generation trend for each turbine (in MWh, Right). The industrial data is scaled slightly to maintain confidentiality. Figure 10. Wake and terrain effects : comparison of power production per turbine for same turbine layout for a. simulated flat terrain, b. simulated complex terrain and c. observed (i.e. measured at field) complex terrain. The industrial observed (measured) data is scaled to maintain confidentiality but the scaling does not influence trends and conclusions.
The comparison reveals that the simulated power production for all turbines for the flatterrain case is substantially lower than that for the realistic complex terrain situations (both simulated and observed complex terrain results). It is important to note that all the simulations are conducted for a single dominant wind-direction at neutral atmospheric condition. The reason for lower power production for flatter terrain is due to two dominant reasons as seen in Figure  11, i.e. (1) due to lack of an ascending terrain induced wind-speed up and (2) due to an extended wake region, which adversely affects power production of those downstream turbines that have their hub at similar altitudes as the upstream turbines. The shorter wake-region in complex terrain can be explained from our observation on wake-terrain interactions, i.e. the terrain induces local turbulence that results in an increased dissipation of the far wake structure. The faster wake-decay in terrain conditions leads to better power production in some downstream turbines. For example, turbine 15 is directly downstream of turbine 14 and their hubs are almost at similar height (around 370-371 m altitude, see figure 11). As a result of slow wake dissipation in flatter terrain case (as seen in figure 9) , wakes from turbine 14 are more stronger when they reach this turbine 15. This results in a larger dip in power production at turbine 15 than turbine 14 for flatter terrain case than for both the complex terrain cases (simulated and observed results) as seen in figure 10). Even though the power produced in turbines 14-15 and all other turbines are higher in complex terrain case due to wind speed-up, the dips in power production from turbine14-15 is more in flatter case , which could be attributed to the wake effects. The comparison also reveals that the simulated complex terrain power generation is higher than the actual annual observed power at the farm-site . The scaled observed averaged annual power generation data for each turbine (in MWh) for the Bessaker wind farm is shown in figure 9 and in figure 10. The discrepancy between the observed data and simulation is expected due to the following reasons , a) the dominant reason being that the measurement (field data observation) involves influence of all wind-directions and all atmospheric conditions over the year. The observed power data accounts for non-productive phases involving stable atmospheric stratification and lower wind-speeds at different directions as compared to the simulated case (which considered the most productive dominant wind-direction with neutral condition), and b) due to scaling of observed data to maintain confidentiality. Since, it was not possible to segregate the observed annual power data for the specific simulation case, so a whole year's observed data was considered for comparison with the simulated flatter terrain and the simulated complex terrain case. Both the studies conducted in this work, i.e. one involving single wake-hunt hill interaction study and the Bessaker wind-farm wake-terrain interaction study, implies that wake-region can have a faster recovery by devising a means to obstruct and annihilate the vortices in wake (tipvortex or vortical structures). The obstruction need not be in direct flow path of wake, but close enough to alter the turbulence levels and induce shear within flow regions (to create secondary vortices). Perhaps, the obstruction can be mounted on the turbine support set-up. This will help to clustering together more turbines in close proximity leading to maximum yields and leading to a smaller footprint on the environment. This work is a first step towards understanding the influence of terrain on wind farm wakes, especially in cases where within a short distance (say 10D) itself, the nearby terrain can rise in altitude and be high enough to influence the wake dynamics. In such a terrain, the placement of turbines can then be closer as wake may dissipate faster.

Conclusion
The current work involves two main studies : (1) Single vortex terrain interaction study, which involve interactions of a single vortex with different simplified terrains to understand the fundamental vortex dissipation mechanism, and, (2) Multi-turbine Bessaker wind-farm waketerrain interaction study, which involves understanding influence of wake and terrain on power production of wind farm using LES simulations and observed field data. The single vortex-terrain interaction study reveals a much faster decay of vortex for the hill cases as compared to a flatter-terrain case. Analysis of vorticity, turbulent kinetic energy and turbulent stress data helps to attribute the faster vortex decay to : (a) the formation of secondary vortices (SV) due to high-shear regions in flow and interaction of these SV with the primary vortex and (b) due to enhanced vorticity diffusion.
The multi-turbine Bessaker wind-farm wake-terrain interaction study reveals that terrain impacts wake dissipation and changes power production.The flatter-terrain has less power production than the complex terrain due to (1) lack of an ascending terrain induced windspeed up and (2) due to an delayed turbine wake dissipation which adversely affects power production of downstream turbines.
Thus, both the studies conducted in this work imply that wake-region can have a faster recovery by devising a means to obstruct and annihilate the vortices in wake (tip-vortex or vortical structures). The obstruction need not be in direct flow path of wake, but close enough to alter the turbulence levels and induce shear within flow regions (to create secondary vortices). This analysis is expected to serve as an initial first step towards forming guidelines for interturbine locations in complex terrain (depending upon the local terrain) and possibility of terrain modification to cause faster wake decay. Future work should involve studying such interactions using practical artificial obstacles in actual wake conditions.