Numerical Study for HAWT Wake shape with different Angles of Attack

Increasing world demand for renewable energy resources as wind energy was one of the goals behind research optimization of energy production from wind farms. Wake is one of the important phenomena in this field. This paper focuses on understanding the effect of angle of attack (α) on wake characteristics behind single horizontal axis wind turbines (HAWT). This was done by design three rotors different from each other in value of α used in the rotor design process. Values of α were (4.8 ̊,9.5 ̊,19 ̊). The numerical simulations were conducted using Ansys Workbench 19Fluent code; the used turbulence model was (k-ω SST). The results showed that best value for extracted wind energy was at α=19 ̊, spread distance of wake behind single HAWT was inversely proportional with α value. Highest turbulence intensity level was in small values of α =4.8 ֯ and 9.5 ֯, which explain presence of meandering phenomenon in these two cases on contrast of α=19 ֯, while symmetry in wake shape was more obvious in α=19 ֯ than other cases due to the moderate turbulence intensity which was found in this case, which indicates that lowest aerodynamics loads on rotor components was in α=19 ̊. So best distribution for turbines in farm could be achieved at use high values of α in comparison to low values.


INTRODUCTION
Request for electricity increases day after day, especially with a huge increase in population and infrastructures. Challenges and risks of extraction, manufacturing and production pollutants make usual energy sources like oil, naturalist gas, and nuclear fuel to produce electricity more complicated than technology advancement in clean renewable energy sources as such. During last years, many experiments and researches in wind energy applications were about the development of wind turbine components as blades, optimization of wind farms design, effects of terrains and topography on wind flow in farms, use numerical simulation software in wind energy applications calculations, losses inherent to extract and generate wind energy like losses from wake phenomenon. Wake phenomenon emerges from upstream wind turbines, and turbine towers shadow in wind ranches, these shapes of wakes overlap with each other. Its impact becomes great on downstream wind turbine rows, which are not limited to power losses but include fatigue damages on rotor blades, tower, and drive train due to aerodynamically loads, which is doubled due to high turbulence levels existing in the wake area. These factors increase the cost of turbine maintenance. Studies were done in large onshore wind ranches indicating that the average extracted power loss because of wakes is approximately 10% -20% of yearly power production (Cole, 2012). So it is very important in the wind farm design process to accurately calculate wake shapes behind wind turbines. There are two interesting critical problems in the wake region: the high level of turbulence intensity and flow speed shortfall, the significant losses in power because of wake and distribution of aerodynamically loads are computed in two regions inside wake space, referred to them: near wake region and far wake region (Adaramola and Krogstad, 2011). The near wake region is characterized by a varying shortfall in pressure and axial velocity attributed to axial thrust (Gómez-Elvira, et al., 2005). Flow at this region is significantly influenced by rotor geometry; the main reason presence of blade tip vortices and sharp changes in pressure and axial flow speed values which cause wake expansion, rotor geometry effects on far wake region are limited on increase turbulence intensity and decrease in wind velocity which is closer to the Gaussian shape (Barthelmie, et al., 2006). The region far away behind the turbine is influenced by atmospheric turbulence and induced turbulence by rotor presence. All these phenomena, characteristics, and features of wake relate to values of aerodynamics parameters that are used in HAWT rotor blades design. The angle of attack one of these parameters, so this study focused on getting a deep understanding of the effect of changing the angle of attack through the design process on the Journal of Engineering Volume 27 June 2021 Number 6 3 behavior of the generated wake and drawing the map of wake achieved for a single wind turbine, which impacts more on the distribution of the wind turbines in the farm.

NUMERICAL ANALYSIS
Computational Fluid Dynamics is an essential numerical method for simulating fluid inflow shapes by resolving equations of fluid motion at a discrete set of points (Al-Zaidee and Kasim, 2018). The used CFD package in this study was fluent code, and this code utilizes the finite volume method (FVM) to resolve fluid flow's partial differential equations. CFD provides various physical models as an in-viscid or viscous, compressible or incompressible, turbulent or laminar inflow (Amori and Adeeb,2016).

Governing Equations
The incompressible Navier stokes equations are applied to describe wake behind the wind turbine, but with excluding for region of the blade tip, two governing equations in Cartesian coordinates and Einstein notation could be used where velocities of atmospheric inflow at the front and behind of wind turbine in the range 4-25 m/s. (Celik,1999) Where (x) is position vector, (μ) viscosity, (u) velocity, (ρ) fluid density, and (p) indicate pressure.
(Sji)strain rate tensor,(t) time, j, and i are directional components, f i external body forces, and defined as:

(RANS) Equations
RANS equations are utilized to deal with turbulent flows; these time-averaged equations for fluid flow motion can be applied with approaches based on flow turbulence properties to obtain an approximate average solution for the RANS equations. These equations are written as: u i ′u j ′,u i u j indicate the fluctuating, instantaneous velocity terms (Celik,1999).

Design Computational Domain
Design of computational domain was done through the Design Modeler, CAD interface in Ansys workbench, the cuboid domain was designed with dimensions (4D, 18D, 4D) as Cartesians (x, y, z) of computational domain represent vertical, streamwise, and spanwise directions, respectively, as depicted in Fig.1. the dimensions of the domain were determined enough for the produced wake and flow rotation. The turbine rotor was located at 2D from the inlet, where D is the turbine rotor diameter. Full HAWT's rotor was enclosed in a cylindrical shape with slightly larger dimensions than rotor diameter and hub thick to avoid any overlap among interface zones, making problems in the Mesh process and errors in numerical simulation.

Grid Generation
The meshing process for the stationary cuboid computational domain involved rotating cylindrical domains around the wind turbine rotor to facilitate the sliding mesh calculation, the meshing process was done through ANSYS Fluent code, and accurate results for the suitable grid (CFD) solutions was generated. An unstructured grid approach was carried out for mesh. This approach is regularly applied infinite volume calculations. Besides its capability to solve complex geometric topologies very well, the non-conformal high-resolution unstructured mesh was made by applying advanced size options, a mixture grid of unstructured wedge, and an unstructured tetrahedral was generated. A precise mesh grid was done for the stationary part domain, as shown in Fig.2. To get accurate resolves for the downstream wake, for best acquiring of inflow behavior in wake area and generated forces near and on the rotor, the high-resolution mesh was done for rotor as shown in Fig.3 and Fig.4. Also, the high-resolution mesh was made for tips of the rotor blades to obtain the best simulation for adjacent boundary layer and emerged vortices from the blades tips, the flow around the blades was captured using inflation layer at value 20 and 1.2 for growth rate option. The first cell was positioned near blades y+ value was selected less than one for all inflow areas surrounding the wall to match the utilized turbulence model, which requires excellent mesh settings close to the wall (Lanzafame, et al., 2014). Table 1 shows parameters of unstructured grid mesh.

Setup and Solution
Three dimensions, double-precision, parallel option, and the solver settings were selected as an incompressible pressure-based solver, absolute value for velocity formulation option, and transient for time. Settings in the model option were determined as viscous k-ω shear-stress transport (SST) closure turbulence model. Solutions of (RANS) were done through the coupled scheme and green gausses node-based option for spatial discretization, with second order option for momentum equation, pressure, turbulence kinetic energy (k), specific dissipation rate (ω), and transient formulation. The used turbulence model had proven its ability to present worthy results in separated areas with a very good guess for the transition process from laminar to turbulent flow (Gómez- Elvira, et al., 2005). Transport equations for k-ω SST model are: Where Ḡ k is the term of production for the kinetic energy of turbulence from average speed gradients, G ω is term of production of ω. Γω and Г k are the effective diffusivity of ω and k, on respectively.Y ω and Y k represent the dissipation ω and k due to turbulence. S k and S ω are userdefined source terms. D ω is defined as cross-diffusion term, where it is the last blending term of the diffusion between two relevant k-ω and k-ε models (Cole,2012). Two types of measurements were performed at different downstream locations behind HAWT position; mean local velocity and turbulence intensity. Measurements locations were 2D, 3D, 5D, 7D, 9D, 12D, where D here is rotor diameter. To save computational time and consideration, the wind turbine tower wasn't involved in turbine design and simulation using computer resources. Specification of used computer was HP Z-book 15 model, workstation core™ i7 processor, RAM 32 GB, video card NVidia 2GB Quadro K2100M,500GB Hard disk.

Boundary Conditions
The velocity inlet option was considered for wind entrance upstream of the HAWT. The outflow option was set for wind exit downstream of the turbine, upper and lateral sides of the computational domain, while the wall option was set for the ground side. The time step size was set at 10 -3 . The inlet turbulence intensity and average wind speed profiles were introduced as 5 % and 5.6 m/s. Each simulation case was done with a 1000-time step.

RESULTS AND DISCUSSION
The horizontal line of wake measurements in the height of the hub was determined for the axial measurement stations (2D, 3D, 5D, 7D, 9D, 12D) downstream of the HAWT.

Wake Characteristics behind with Different Angles of Attack
Comparative analyses were made depending on normalizing the downstream profile distance at the Y-axis, vertical profile distance at the X-axis, and lateral profile distance Z-axis with respect to the diameter of the wind turbine rotor. Fig.5 represents the lateral profile of wake (YZ plane) at the hub center. It shows that length of the near wake is 5.7D at case α = 4.8 ֯ and 5.3D and 3.2D at case α =9.5 ֯ and case α=19˚ on correspondingly, and linear relation was found between length of near wake and value of α: Also, the maximum width for the wake plane was adjacent to HAWT in all cases, and the amount of maximum width of the wake was noted at case α=4.8 ֯ with value 1.38D. The smallest width presents at α=19֯ with value 1.2D, which means presence inversely proportional relation for wake shape width with the value of α in lateral profile. Also, a high value for velocity deficit surrounded by severe velocity gradient was observed in this profile. Still, this deficit gradually becomes less sharp further downstream and then finishes due to momentum transportation from the freestream into the wake area. The phenomenon of meander waking, which plays a vital role in increase the aerodynamics loads on downstream turbines in farm and disturbs the stability of wake profile in the far area of the flow field, was clearly observed at two cases α=4.8 ֯ and α=9.5 ֯ in comparison to case α=19֯ . The strength of this phenomenon is related to turbulence intensity as indicated in previous researches, so this explains the reason behind the simple meanderings in case α=19֯ , that turbulence intensity in this case, as clear in Fig.6 where the lowest value among three cases of α. Semi-symmetry in wake shape in three cases of α was found and extended until 3.8D, 3.3D, 3.2D in cases α=4.8 , α=9.5 ֯ and α=19 ֯ on respectively. Maximum local velocity deficit was presented in cases α=4.8˚ and α=9.5 ֯ with uniform style, while the lowest value was observed in case α=19˚ but with non-uniform style. Fig. 6 represents the vertical profile of wake (XY plane) at the hub center. It shows that the length of wake path is about 5.44D in case of α = 4.8˚ and case of α =9.5˚, while in case of α =19 ֯ is 3.8D. The semi -symmetry in near wake shape in three cases was found. About meandering in wake movement, the most stable wake was in case α=19˚ in comparison to the two cases, where the meandering phenomenon was clear, especially in case α=4.8˚. Regarding wake's width value, a small difference was observed in the three cases with values 1.6D,1.8D,1.75D in α = 4.8֯ , α = 9.5 ֯ and α = 19 ֯ on respectively.Maximum velocity deficit was demonstrated at α = 4.8֯ adjacent to the upper and bottom sides of the rotor at the tips of the rotor. This deficit extended until axial distance 4.2D, same style was showed at case α = 9.5֯ but vanished at 4D. In contrast, in case α = 19 ֯, the lowest value for velocity deficit was observed among the three cases. The wake structures behind a single wind turbine were analyzed using the turbulence intensity contours at selected streamwise stations (Y/D = 2, 3, 5, 7, 9, 12). At station 2D, Fig.7 indicates that at case α=19 ֯, the high turbulence intensity value was (9%) with spreading to the upper side more than the bottom side of the rotor. In comparison, at cases α=4.8֯ and 9.5֯ highest turbulence intensity (11.7 % -12 %) were observed as distinctive ring or frame behind rotor blade tips which due to highly turbulent shear layer shaped by vortices of blades tips, there were three un-equally cores of this value observed at case α=4.8֯ and two cores at case α=9.5֯ . More downstream, at station 3D behind the wind turbine, the turbulence intensity level was slightly less than the station 2D. Still, its diffusion was wider to include the rotor swept area as observed in three cases in Fig.8. At station 5D, as indicated in Fig.9, the wake turbulence intensity level in case α=4.8֯ was characterized by a core of high turbulence level (10% ) at the upper part of rotor swept area, then started in decrease and spreads more to topside 1.1D and left side 0.85D of rotor swept area, while in case α=9.5֯ turbulence intensity level was less than the case of α=4.8֯ , in range (7% -9%) and extended to 1.2D at the upper side and 1.35D at the left side of rotor swept area, while in case α =19֯ lowest value of turbulence intensity was observed in range (4%-6 %) in comparison with other cases. Still, its vertical spread to the upper side was the biggest 2.4D, with spanwise extend 1.1D to the right side. For station area 7D at case α=4.8 as indicated in Fig.10, turbulence intensity profile continued in its spread to upper and left parts from field flow, that turbulence intensity was stronger at left part than the upper part, the difference in case α=9.5֯ was spread to right and left sides to add to the upper side, but spread area was less in the upper side and wider on spanwise sides. Still, turbulence intensity value was in same range of case α=4.8֯ , in case α=19֯ clear decrease in level of turbulence intensity was observed in comparison to other cases of α. The spread area was narrower in the spanwise direction. At station 9D behind wind turbine Fig.11 turbulent intensity value decreased to a moderate level at case α=4.8֯ with diffusion to upper side to more than 2.25D in x-axis with semi-symmetry, while at case α=9.5֯ moderate level of intensity of turbulence was revealed but with meandered movement in turbulence profile, and asymmetry shape moved to upper side till 1.85D, at case α=19.֯ Also, a moderate level of turbulence intensity was shown but with a symmetrical profile around the x-axis extended to the upper side until 1.44D.In 12D downstream station as indicated Fig.12 symmetry in turbulence profile around x-axis at two cases, α=4.8֯ and α=9.5֯ was clear, and distinguished by two cores of moderate level of turbulence intensity on sides of rotor swept area or Journal of Engineering Volume 27 June 2021 Number 6 10 blades' tips tracks, which is extended in upper side to 1.5D Also, symmetry was recognized in turbulence profile around x-axis at case α=19֯ but with less turbulence intensity level expanded in the vertical direction to 2.1D, which mean that turbulence profile in this station was largest in case α=19֯ in comparison to the two other cases. Extracted power in each of three cases, was 0.566 watt in case α=4.8֯ , 0.337 watts, 1.58 watt in cases α=9.5˚ and α=19֯ on respectively, which indicate that the highest extracted power among the three cases was at case α=19֯ .

CONCLUSIONS
Three different angles of attack were used in the design processes of three rotors. By numerical solution, the investigation was done for the relation between the wake characteristics and α value, and possibly applies it to design the best wind turbine distribution map in the wind farm. Numerical simulation results revealed that the spread of near wake in lateral and vertical profiles of the flow field was inversely proportional with α value, which is attributed to rapid decay for the shear layer of vortices structure, which be comprised from primary vortices of rotor blade tips and secondary vortices of the hub. The meandering phenomenon in wake path at lateral and vertical profiles was more clear in cases α=4.8֯ and α=9.5֯ than case α=19֯ , which be attributed to high turbulence intensity in these cases of α. For maximum width of wake shape in lateral profile was same in cases α=4.8֯ , α=9.5֯ , and the lowest value was in the third case, while in vertical profile, the values of maximum width for wake were close in three cases, about the symmetry of wake shape in lateral and vertical profiles was noted that symmetry phenomenon was the case in α=4.8֯ for a longer distance than other two cases, and symmetry in vertical profile extended more than symmetry in lateral profile for the three studied cases of α. About turbulence intensity in (ZX) plane, the highest levels were observed in cases α=4.8֯ and 9.5֯ , with further movement downstream, the turbulence intensity decreased. Its shape expanded more on sides of rotor swept area accompanied by the movement to the upper side of the flow field, also symmetry in turbulence intensity shape in this plane was found at some of the measurement stations in three cases of α. For extracted power, the highest value was in case α=19֯ . All findings in this study indicated to importance of design HAWT rotor blades by the high value of α, which provide best-extracted power from wind and fewer aerodynamics loads on turbines in the farm, in addition to possible use of small or limited area to design wind farm with good production