Numerical simulation of flow over an airfoil for small wind turbines using the γ-Reθ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma - {\text{Re}}_{\theta }$$\end{document} model

The mechanism of the laminar separation bubble and the laminar-turbulent transition over the airfoil UBD5494 is simulated in ANSYS-FLUENT using the transition gamma-Reθ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{gamma}} - {\text{Re}}_{\theta }$$\end{document} model at Reynolds number 6 × 104, 1 × 105, 1.5 × 105, 2 × 105 and 3 × 105. Modified constants of the Reynolds momentum thickness are incorporated in the model. The aerodynamic performance of the airfoil is also examined against the flow behaviour. Simulation results show that with the increase in angle of attack, laminar separation bubble moves towards the leading edge and at the same time contracts in size. It starts to expand after reaching the foremost point of the leading edge and then bursts, resulting in flow turbulence and stall. With decreasing Re, the size of the laminar separation bubble is found to be increasing and its progress towards the leading edge is noticed to be slower. The numerical results also indicated that UBD5494 airfoil has enhanced lift-to-drag ratio and desirable stall characteristics which are distinctively advantageous for the better performance of small wind turbine rotors.

Specific turbulence dissipation rate (m 2 s -1 ) k Turbulent kinetic energy (J kg -1 ) S Strain rate (s -1 ) y Distance from to nearest wall (m) y?
Distance in wall coordinates

Introduction
With the increasing world's energy requirement and growing environmental concerns, renewable energy resources like wind, solar, biomass, etc. are expected to significantly supplement the expensive and depleting fossil fuels. Among these, wind power with its cumulative installed capacity of nearly 31 9 10 5 MW by 2013 has alone contributed about 2.5 % of worldwide electricity demand [1]. It is expected to further grow and reach up to 8-12 % by 2020 [1]. In contrast to large wind power sector, small wind industry is also growing at a fast rate. With the growth rate of nearly 35 % annually, total installations of small wind turbines (SWTs) in 2015 are projected to be approximately 400 MW [3]. These installations include both the Horizontal and Vertical axis type wind turbines (HAWTs and VAWTs). From 2015 to 2020, with 1000 MW of newly installed capacity added annually, a steady growth rate of 20 % is forecasted, leading to a cumulative installed capacity of 5 GW by 2020 [3].
SWTs are defined as systems with rotor swept area not more than 200 m 2 with an equivalent power of about 50 kW [2]. Applications of such small wind systems include, but are not limited to the sectors of water pumping, telecommunication power supply, irrigation, homes and small industries [4]. Installations of such small systems are often based on the places where the power is required and not necessarily on the strength of the wind [5,6]. Additionally, for simplicity reasons, they are designed to selfstart and work as stand alone or in connection with microgrids [5,6].
Moreover, due to the small blade size and the low wind speed working conditions, airfoil device employed along the small HAWT blades operates at low Reynolds number (Re) from hub to tip [6]. The complex nature of flow about the blades in low Re warrants for a careful and 'clever' selection of the airfoil profiles in the SWT design process [7]. Typical airfoils designed for high Re, such as NACA airfoil series, are reported to underperform under such low Re conditions and consequently degrade the wind turbine performance [7,8]. Alongside, airfoils designed for low Re specifically for small HAWTs are found to be limited in the literature [7][8][9][10][11][12][13][14][15][16][17].
While developing small HAWTs for specific applications, the initial aerodynamic performance of an airfoil is generally investigated experimentally in a wind tunnel. At low Re, such wind tunnel measurements are reported to be challenging due to the requirement of higher level of accuracy in equipments for correct modelling of the flow around the airfoils [18]. These include the modelling of the laminar separation bubble (LSB) and the flow behaviour over airfoil's surface. For example, under the same testing conditions, a difference of 50 % in drag coefficient measurement of Wortmann FX63-137 airfoil was reported under three different testing facilities [19]. On the other hand, with modern computational power, it is now possible to simulate transitions and turbulent flows with Reynolds averaged Navier-Stokes (RANS)based computational fluid dynamics (CFD) models with lower risks of inaccuracy [19,20]. Yet, a precise computational study requires proper modelling and interpretation of the transition physics. This paper is formulated by keeping such level of accuracy in the simulation and modelling in mind.
Moreover, this paper is a first computational effort towards the understanding of the flow behaviour about the airfoil UBD5494 at low Re for small HAWT applications, using the Menter's c -Re h model. The model chosen has distinct advantage of associating the transition model with experimental data [21][22][23][24]. Much attention has been paid to the understanding of the LSB and its direct consequences on the airfoils aerodynamic performance. The results from this study will be supportive for the further development of new low Re airfoils.
A brief description of the LSB is given in Sect. 2, which is followed by the features of the chosen airfoil in Sect. 3. Computational setup and transition model is described in detail in Sect. 4. Section 5 provides a discussion on computational results. Finally, conclusions are drawn in Sect. 6.

Laminar separation bubble (LSB)
Reynolds number (Re), which is the ratio of inertia forces over viscous forces, is considered to be low if Re \ 5 9 10 5 [6][7][8][9][10][11][12][13][14][15][16][17][18]25]. Airfoils operating under such low Re, at certain range of angles of attack, a (measured in degree), face separation in the attached boundary layer due to the inability of the flow to overcome the adverse pressure gradients (APG). Under this situation, the transition from laminar to turbulent flow occurs in the free shear layer, followed by reattachment, which forms the so-called LSB, as shown in Fig. 1.
The size of the LSB is defined by the expanse between the point of flow separation and reattachment (Figs. 1, 2) and is generally expressed in terms of the airfoil chord length (x/C) (measured in meter). The airfoil's chord length (x/C) is the distance between the trailing edge of the airfoil and the point on the leading edge where the chord intersects the leading edge. At Re = 1 9 10 5 , longer LSB is generally reported to cover up to x/C = 20-30 % of upper surface of airfoil [25], whereas at higher Re, the LSB can be shorter and harmless [25]. LSB can be interpreted from the hump in the actual pressure distribution plot (solid line) about the airfoil, as shown in the pressure distribution plot in Fig. 2 [8]. Pressure distribution plot represents the pressure at each point about the airfoil, at the given a, and is represented by the dimensionless number, the pressure coefficient (C p ). For inviscid flow (dashed line) about the airfoil, shown in Fig. 2, no LSB will be observed [18].
The mechanism of LSB was first observed by Jones in 1930s [26] and further explored by Gaster in 1960s [27]. They analysed the stability behaviour of the separation bubble experimentally. A novel semi-empirical bubble model was developed by Horton in 1968 [28] which has been broadly employed in such studies. Further contributions in the characterization of LSB were made by McGregor [29], Young and Horton [30] and Woodward [31]. Further, extensive literature on LSB and low Reynolds aerodynamics was documented by Carmichael [20], Lissaman [25] and Tani [32].
The formation of LSB greatly influences the skin friction drag (C f ). The skin friction drag (C f ) can be expressed by the von Kármán integral boundary-layer equation [26] as where n is the boundary-layer co-ordinate, u e is the boundary-layer edge velocity. The term H is shape factor, defined as the ratio between the boundary-layer momentum thickness h ð Þ and the boundary-layer displacement thickness d Ã ð Þ. The boundary-layer displacement thickness h ð Þ is defined as the distance perpendicular to the boundarylayer surface to which the boundary surface has to be exiled outwards to balance the decrease in the discharge owing to the configuration of the boundary layer [34]. Whereas momentum thickness d Ã ð Þ is the distance measured perpendicular to the boundary surface, to which the boundary has to be displaced for the reduction in momentum of flowing fluid as it is reduced by the existence of the boundary bubble [34].
The skin friction inside the bubble is nearly zero [8] and hence the integral boundary condition becomes where qu 2 e h represents the total momentum defect. Considering the boundary-layer edge velocity, the shape factor and that the total momentum defect as average quantities, the above equation can be re-written as It can be seen that the increase in drag D qu 2 e h À Á due to LSB is proportional to the product of average mass defect qu e d Ã and edge velocity jump Du e . Here it is to emphasize that  the size and position of the LSB are the functions of the airfoil shape, angle of attack (a), Re and environmental interruptions [28,32].

UBD5494 airfoil
Airfoil UBD5494 was specifically designed to enhance the performance in terms of lift-to-drag (L/D) ratio. A 3D view of the airfoil is shown in Fig. 3. It was designed to be employed over the entire blade of a small HAWT with rated power B1 kW. However, for small HAWTs with rated power greater than 1 kW, UBD5494 is recommended for the tip region only. The direct design method was used to extract the final shape of UBD5494 from an existing airfoil Eppler-62. X-foil, a post-design viscous/inviscid analysis tool, was used for the designing process. As detailed by the name, UBD5494 owns maximum thickness of 4.81 % at x/C = 22.8 % and camber of 5.4 % at x/C = 49.1 %. For further design specification see [7]. UBD5494 has been reported to produce higher L/D ratio compared to the other existing airfoils while operating at low Re and thus enhancing the performance of small HAWTs [7].

Computational method
Computational domain and solver setup A structured grid 2D C-H topology quadrilateral mesh was generated around UBD5494 airfoil in the mesh generating tool ICEM-CFD. Quadrilateral-type cell was chosen because it can provide high-quality solution with less number of cells compared to the triangle mesh [24]. Moreover, in order to attain a fully developed and expanded flow, the length of the computational domain was made 40 times that of the chord length of the airfoil, whereas the width is kept 30 times. To ensure the computed aerodynamic results are independent of the grid size, the density of grid was increased until negligible difference in solution is attained towards convergence. Such methodology of optimum grid selection is adapted from [35,36]. At the outset, the coarse gird named as Grid 1 containing 50,365 cells was made. Number of cells was increased to 70,827 in Grid 2, whereas Grid 3 and Grid 4 contained 108,562 and 130,204 cells, respectively. For all the grids, the mesh density in the attached boundary layer was increased to capture the transition, flow separation and most importantly the predicted separation bubble. The mesh density was kept progressively coarser in the far-field area where the flow gradients approaches zero. Dense grids, with increased number of cells, were also placed near the leading and trailing edge because of the steepest gradients. The transition in mesh size was kept as smooth as possible for numerical accuracy. The boundary conditions for the inlet were taken as velocity inlet and the outlet conditions were defined as the pressure outlet.
The solver was set for steady state. Desired angle of attack, a, was attained by rotating the mesh in order to capture the movement of possible separation bubble and the stall angle. Velocity at inlet was specified to achieve the desired Re. To solve the coupled problem between pressure in momentum equations and velocity components, semiimplicit method for pressure-linked equations (SIMPLE) algorithm [24] was employed, and second-order upwind spatial discretization was set in calculation. The spatial gradient was selected as the least squares cell based. Air pressure was taken as standard. Convergence criteria residual target values were set to 10 -6 . For precise simulation of the boundary-layer flows and the coupled lift and drag forces of airfoils, y?, which is a non-dimensional distance from the wall to the first node of the mesh was satisfied with the value of y? \1 (Fig. 3).
As shown in Fig. 4a and b, difference in the aerodynamic lift and drag coefficients for a = 7°becomes negligible from Grid 2 onwards. Thus, the refined Grid 3, with 108,562 cells, was adopted for numerical experimentation. Image of Grid 3 is shown in Fig. 5. Final Numerical analysis on the airfoil was performed using ANSYS-Fulent at Re = 6 9 10 4 , 1 9 10 5 , 1.5 9 10 5 , 2 9 10 5 and 3 9 10 5 with the above-mentioned solver settings.

Transport equations
In the present work, transition model c À Re h is used to incorporate turbulence in the flow. It is a four equation turbulence model that combines Shear Stress Transport komega turbulence model (SST k-x) transport equations with two additional transport equations, one for intermittency (c) and second for Transition Reynolds number (Re ht ). Here, intermittency term is employed to activate the production term of the turbulent kinetic energy (TKE), downstream of the transition point in the boundary layer, whereas the Transition Reynolds number term captures the non-local effect of the turbulence intensity [21]. This model is reported to have a distinct advantage of associating transition modelling with experimental data [21][22][23][24]. Moreover, according to [21][22][23][24], the transport equation for intermittency, c; read as where l is the molecular viscosity (Pa.s) and l t is the eddy viscosity (Pa.s). The transition sources, P c1 and E c1 ; are defined as And the destruction sources, P c2 and E c2 ; are defined as where S represents the strain rate (s -1 ) and F lenght is the empirical correlation to control the length of the transition region. The term X represents the vorticity magnitude. The terms c e1 ; c a1 ; c e2 ; c a ; r Ç are constants in the intermittency equation, with values c e1 ¼ 1:0; c a1 ¼ 2:0; c e2 ¼ 50; c a2 ¼ 0:06; c a ¼ 0:5; r Ç ¼ 1:0: The F onset determines when the production of c; the intermittency; is activated. The functions to control transition onset are given by Re hc is the critical Reynolds number where intermittency begins to increase in the boundary layer. This happens upstream of the transition Reynolds numberRe ht and the difference between the two must be obtained from an empirical correlation [21][22][23]. Both the Re hc correlations are functions ofRe ht : The term Re T is the viscosity ratio and Re V is the strain rate (vorticity) Reynolds number, whereas k is the TKE (J kg -1 ) and x represents the specific turbulence dissipation rate (m 2 s -1 ).

Separation-induced transition modification
The tailored separation-induced transition is written as Once the viscosity ratio is large enough to force reattachment, F reattach disables the modification. F ht is the blending function used to turn off the source term in the boundary layer: F wake ensures that the blending function are not active downstream of an airfoil. The transport equation for transition momentum thickness number,Re ht ; is expressed as The source term, P ht ; is defined as Here, The parameter t is determined using dimensional analysis. The values of the constants c ht and r ht in the Reynolds momentum thickness are given in [21][22][23][24] as c ht ¼ 0:03; r ht ¼ 2:0; where term c ht controls the magnitude of source terms and r ht controls the diffusion coefficient. These constants values in the original model were derived on flat plate transition experimentation [23]. However, in order to obtain physically realistic solution, modified constants based on low Re simulations [33] of circular arc airfoils are adapted in the present work, given as c ht ¼ 0:02; r ht ¼ 3:0:

Coupling transition model and SST transport equations
Finally, the transition model interacts with SST transition model as where Y k and G k are the original terms of SST model for destruction and production, respectively. The production term in the x equation was unchanged.

Results and discussion
The performance of UBD5494 airfoil and the flow behaviour over its surface acquired from ANSYS-FLUENT are shown in Fig. 6 through Fig. 12. The variations in C l and C d with a are shown in Fig. 6. It can be seen that for all Re, C l increases with alpha up to a certain angle of attack and then it drops. The initial linearity represents the attached flow over most of the airfoils surface and the drop represents the stall angle at which the flow is mostly detached from the surface. Between these two scenarios, the C l values are observed to be staying quite constant which demonstrates the soft stall behaviour. Such manner of stall is understood to be mainly because of the lower camber and/or moderate thickness along with small round leading edge of the airfoil. For all the Re examined, the airfoil remains in soft stall from a = 7°to 11°, followed by the stall angle at a [ 12°. For small HAWT blades, such stall characteristic is desirable as it does not increase the mechanical loading on the components [6].
As expected, the L/D ratio increases with the increase in Re. For Re = 6 9 10 4 , 1 9 10 5 and 1.5 9 10 5 , the airfoil could produce a maximum L/D ratio of 40, 60 and 75 at

Angle of aƩack, degree
Re= 60 000 Re=100 000 Re=150 000 Re= 200 000 Re=300 000 Coefficient of drag, Cd Coefficient of liŌ, Cl Fig. 6 Coefficient of lift and drag as a function of angle of attack a = 6°, respectively. For higher Re, that is 2 9 10 5 and 3 9 10 5 , the maximum L/D values of 85 and 108 were observed at a = 4°. These maximum L/D values obtained at a = 4°and 6°will be significant in early capture of maximum energy from small HAWTs.
Moreover, Fig. 8 shows the comparison of the distribution of the pressure coefficient at different Re while a is fixed at 7°. As expected, the suction pressure peak increases with increase in Re. A long hump at Re = 6 9 10 4 indicates a longer separation bubble. The point initiating the constant pressure represents the point of separation of laminar boundary layer from the surface, and pressure values keep constant till the transition occurs. The increase in pressure which can be seen after transition to the turbulent flow is mainly due to the turbulent entertainment which energizes the shear layer [28]. The pressure trend between transition and reattachment points is seen to be almost linear and is expected to be accompanied by the growth of the shear layer [8]. With the increase in Re, the LSB moves towards the leading edge of the airfoil. Along with this movement, LSB is also found to contract in size. This is clearly revealed from the extent of separation and reattachment points in Fig. 9.
Moreover, within the LSB, the laminar portion is seen to decrease and the turbulent portion increases in size with the increase in Re. This is further evident from the velocity contours shown in Fig. 9a-e, captured from ANSYS. In the velocity vector plots with constant a = 6°a nd varying Re, the reattachment of the turbulent flow, along with the shortening and movement of LSB can be clearly observed. As the LSB gets shorter, most of its area is covered by reversed flow vortex. The shortening of LSB with increasing Re is mainly caused by the energized flow that overcomes the APGs and forces the turbulent flow to reattach to the airfoil surface. In order to further investigate LSB behaviour, their extents are discerned from the computed skin friction distributions from ANSYS-FLU-ENT and plotted against angle of attack as shown in Fig. 10.
At Re = 6 9 10 4 and a B 4°, the laminar boundary layer remains attached over most of the airfoils suction surface area and leaves the airfoil surface smoothly. As the angle of attack increases, the separation point moves further upstream followed by transition with no reattachment. The calculations show that the bubble has developed into the wake under this condition. From Fig. 10a, the LSB has started to occur at a = 6°, where the laminar boundary layer separates at x/C = 0.5 and reattaches ahead of the trailing edge at x/C = 0.85. At this point, the length LSB is calculated to cover 35 % of the airfoils surface. In addition, it is also noticed that the span of the laminar potion of the bubble (from separation point till the transition point) is longer compared to that of the turbulent region (from the transition point till the reattachment point)-shown in Fig. 9a. The computed TKE past transition for this point is visualized in Fig. 11a, whereas the velocity contour, which visualizes the flow velocity over the airfoil, at this point, is shown in Fig. 11b.
Further from Fig. 10a, at a = 7°, the scenario of separation, transition and reattachment occurs at x/C = 0. 2, 0.3 and 0.42, respectively, whereas at a = 8°and 10°, the bubble further shifts towards the leading edge. Beyond a = 10°, the LSB has started to expand, with decrease in laminar portion and increase in turbulent portion. At a = 12°, the LSB is at the verge of bursting. At a [ 12°, the bubble bursts and the turbulent flow does not reattach to the airfoil surface bringing the airfoil to stall.
In Fig. 10b, flow separation along the airfoil's chord length is plotted as a function of a at Re = 1 9 10 5 . The flow pattern is found similar to Re = 6 9 10 4 and the LSB

Pressure coefficient, Cp
Chord lenght, X/c Re=60 000 Re=100 000 Re=150 000 Re=200 000 Re=300 000 Fig. 8 Pressure distribution plots for UBD5494 airfoil at various Reynolds number for a = 7°I nt J Energy Environ Eng (2015) 6:419-429 425 Fig. 9 Laminar separation bubble formation over UBD5494 airfoil's suction surface a a = 6°at Re = 6 9 104, b a = 6°at Re = 1 9 10 5 , c a = 6°at Re = 1.5 9 10 5 , d a = 6°at Re = 2 9 10 5 and e a = 6°at Re = 3 9 10 5 has started to develop near the trailing edge from a = 6°. The separation, transition and reattachment happen at x/ C = 0.52, 0.63 and 0.8, respectively. At a = 11°, the turbulent flow reattaches further downstream, forming a longer LSB. With further increase in a, the bubble bursts at Re = 1.5 9 10 5 and 2 9 10 5 , and LSB is found to develop till a = 13° (Fig. 10c, d). However, the overall length of the LSB is found to be reduced with the increase in Re. At a = 10°and Re = 3 9 10 5 , the movement of transition towards leading edge with higher Re can be observed from TKE contour and velocity contours near the leading edge, as shown in Fig. 12a and b. Moreover, for wind turbine application, examining the influence of airfoils surface roughness is also essential. From C l -a alpha plots, drag-polar plots and the LSB movement, it is revealed that maximum lift for UBD5494 is occurring while turbulent flow is covering about up to 90 % of the airfoils suction surface. This postulates that the a Re = 6 9 10 4 , b Re = 1 9 10 5 , c Re = 1.5 9 10 5 , d Re = 2 9 10 5 and e Re = 3 9 10 5 airfoil is either not sensitive to roughness or the roughness effects are diminished due the higher turbulence level.

Conclusion
In this study, the performance of the UBD5494 airfoil under different Reynolds numbers is analysed using c À Re h . Transition model. Simulations were done at Reynolds numbers of 6 9 10 4 , 1 9 10 5 , 1.5 9 10 5 , 2 9 10 5 and 3 9 10 5 . At Reynolds number 6 9 10 4 and angle of attack 6°, the LSB is observed to be longer, which covers up to 35 % of the airfoils mid-chord section. The trend of contraction and extension of the bubble are attributed to the small leading edge nose and stronger APGs with increase in angle of attack. It is found that the distance between the points of laminar separation and transition has major influence on the bubble as a whole. Due to very low thickness of the bubble, its effects on the aerodynamic performance of the airfoil are not prominent. Thus, the UBD5494 airfoil is found to show smooth stall behaviour. The overall analysis clearly indicates that the aerodynamic performance of UBD5494 airfoil makes it a potential candidate for the blades of small HAWTs.