Computational Modelling of Rectangular Sub-Boundary Layer Vortex Generators

: Vortex generators (VGs) are increasingly used in the wind turbine manufacture industry as ﬂow control devices to improve rotor blade aerodynamic performance. Nevertheless, VGs may produce excess residual drag in some applications. The so-called sub-boundary layer VGs can provide an effective ﬂow-separation control with lower drag than the conventional VGs. The main objective of this study is to investigate how well the simulations can reproduce the physics of the ﬂow of the primary vortex generated by rectangular sub-boundary layer VGs mounted on a ﬂat plate with a negligible pressure gradient with an angle of attack of the vane to the oncoming ﬂow of β = 18 ◦ . Three devices with aspect ratio values of 2, 2.5 and 3 are qualitatively and quantitatively compared. To that end, computational simulations have been carried out using the RANS (Reynolds averaged Navier–Stokes) method and at Reynolds number Re = 2600 based on the boundary layer momentum thickness θ at the VG position. The computational results show good agreement with the experimental data provided by the Advanced Aerodynamic Tools of Large Rotors (AVATAR) European project for the development and validation of aerodynamic models. Finally, the results indicate that the highest VG seems to be more suitable for separation control applications.


Introduction
In order to achieve the aim of 100% renewable energy consumption, wind energy, as a leader toward the way of renewable energy, is developing rapidly all over the world.To decrease the levelized cost of energy (LCOE), the size of a single wind turbine has been increased to 10 MW nowadays, and it will increase further in the near future.Large wind turbines and their related wind farms have many challenges in aerodynamics, aero-elasticity and aero-acoustics.Therefore, the introduction of new designs and upgrading methods will be also needed to make this possible [1,2].
According to Aramendia et al. [3,4] and Fernandez-Gamiz et al. [5], either the use of passive controllers (vortex generators (VGs), microtabs, spoilers, etc.) or active controllers (trailing-edge flaps, synthetic jets and air jet vortex generators) can be introduced.One of the best options is the introduction of VGs on the blades to improve the aerodynamic performance needed.The use of this passive flow device helps to reduce the load distribution, and the blade has to undergo and obtain the optimal creation of energy.Øye [6] compared the measured power curves with VGs and without on a 1-MW wind turbine.Although quite rough methods were used for the VG design optimization, the experiment showed that, for a stall-regulated wind turbine, power increased nearly 24% by using VGs through field tests.Furthermore, Sullivan [7] conducted an experiment on a 2.5-MW wind turbine to test the effects of adding VGs on the power conversion performance.An increase of 11% in the annual energy production was found.
In the work of Fernandez-Gamiz et al. [8], Blade Element Momentum BEM-based computations were carried out on the National Renewable Energy Laboratory (NREL) 5-MW baseline wind turbine with and without flow control devices.The results obtained from the clean wind turbine without any device for flow controlling were compared with the ones obtained from the wind turbine equipped with vortex generators and Gurney flaps.A best configuration case is proposed, which has the largest increase of the average power output.In that case, increments on the average power output of 10.4% and 3.5% have been found at two different wind speed realizations VGs were first introduced by Taylor [9], and their principle of operation relies on the increased mixing between the external stream and the boundary layer due to longitudinal vortices produced by the VGs.Fluid particles with high momentum in the streamwise direction mix with the low-momentum viscous flow inside the boundary layer; therefore, the mean streamwise momentum of the fluid particles in the boundary layer is increased.The process provides a continuous source of momentum to counter the natural boundary layer momentum decrease and the growth of its thickness caused by viscous friction and adverse pressure gradients (Doerffer et al. [10] and Steijl et al. [11]).Vortex generators can reduce or eliminate flow separation in moderate adverse pressure gradient environments (Velte et al. [12,13]).Even when separation does occur for cases of a large adverse pressure gradient, the mixing action of trailing vortices will restrict the reversed flow region in the shear layer and help maintain some pressure recovery along the separated flow.Thus, the effects of separation may be minimized or even removed.
Fernandez-Gamiz et al. [14] and Urkiola et al. [15] studied the behavior of a rectangular VG on a flat plane and the streamwise vortices produced by them to investigate how the physics of the wake behind VGs in a negligible streamwise pressure gradient flow can be reproduced in Computational Fluid Dynamics CFD simulations and their accuracy in comparison with experimental observations.In the work carried out by Gao et al. [16] and Baldacchino et al. [17] on a 30% thick DU97-W-300 airfoil (Delft University of Technology, Delft, The Nederlands), the maximum lift coefficient was substantially increased due to the implementation of passive VGs.When the angle of attack increases, both lift and drag coefficients rise up to values higher than the ones reached in steady state conditions.
The concept of micro vortex generators was most probably first introduced by Keuthe [18].In that study, wave-type micro VGs with a height of 27% and 42% of the local boundary layer (BL) thickness were mounted on an airfoil to decrease trailing edge noise by suppressing the formation of a Karman vortex street and by reducing the velocity deficit in the airfoil wake.Since the late 1980s, these devices appeared in the literature under different names such as sub-boundary layer vortex generator (SBVG), presented by Holmes et al. [19], the submerged vortex generator of Lin et al. [20], the low-profile vortex generator Martinez-Filgueira et al. [21] and the micro vortex generator (Lin [22]).As stated in the review of Kenning et al. [23], the potential applications of VGs and SBVGs include control of leading edge separation, shock-induced separation and smooth surface separation.
Vortex generators are applied on wind turbine blades with the major aim to delay or prevent the separation of the flow and to decrease the roughness sensitivity of the blade.They are usually mounted in a spanwise array on the suction side of the blade and have the advantage that they can be added as a post-production fix to blades that do not perform as expected.Therefore, adding VGs is a straightforward solution to improve the rotor performance (Bragg et al. [24]).
The main contribution of this work is the computational study of the primary vortex generated by three different rectangular VGs, where it has been found that the highest vane is the most suitable actuator for flow control applications.The goal of this study is to investigate how well the simulations can reproduce the physics of the flow behind a rectangular VG to characterize the primary vortex generated by the vane mounted on a flat plate with three different device heights.The aspect ratio (AR) of the vane defined as the relation between the height and length of the VG will consequently vary, maintaining the angle of attack constant.The baseline VG has an aspect ratio of AR H = 2.5, and two height variations are also investigated, those corresponding to AR H2 = 2 and AR H1 = 3, as sketched in Figure 1.For this purpose, computational simulations have been carried out using the RANS (Reynolds averaged Navier-Stokes) method and at a Reynolds number Re = 2600.The case consists of a single vortex generator on a flat plate with the angle of attack of the vane to the oncoming flow β = 18 • .The flow over the flat plate without the VG has been previously simulated to calculate the boundary layer thickness at the VG position.

Baseline Experimental Data
The current study is based on Advanced Aerodynamic Tools of Large Rotors (AVATAR) European project, and the parameters used for the validation of the computational results are the following ones:

•
Streamwise velocity profiles at different locations (AVATAR Task 3.1 [25]) Velocity fields at streamwise planes 10H, 25H and 50H (AVATAR Task 3.1 [25]) Fields of the normalized vorticity at different streamwise planes (AVATAR Task 3.2, [26]) Turbulence kinetic energy fields at streamwise planes (AVATAR Task 3.2 [26]) Taking into account the previous configuration and in order to characterize the primary vortex generated by the three different VGs, the following parameters have been proposed:

•
Streamwise velocity profiles at different streamwise locations.

•
Primary vortex vertical and lateral positions to define the vortex trajectory.
According to the experimental data available in the AVATAR project and to reproduce the experiments of the VG on a flat plate performed in that project, the numerical simulations have been carried out at a Reynolds number Re = 2600 based on local BL momentum thickness θ = 2.4 mm and with a free stream velocity of U∞ = 15m/s.A negligible pressure gradient on the plate is assumed.The angle of incidence of the vane to the oncoming flow is 18°.The calculation of the local BL momentum thickness θ was made by the application of Equation (1): where ux is the streamwise velocity component, U∞ the free stream velocity and y the vertical coordinate normal to the wall.In order to validate the numerical model of the VG, the results of the numerical simulations are compared to the experimental data of the AVATAR Project [25,26].The experiments were carried out in the Boundary Layer Wind Tunnel of the Delft University of Technology.The tunnel can attain a maximum speed of 38 m/s in the wide-walled test section of 1.5 × 0.25 m 2 .The large separation of the side walls (1.5 m) minimizes end effects on the flow region of interest along the centerline zone.An adjustable back wall allows adjustment of the pressure gradient, ensuring a truly null pressure

Baseline Experimental Data
The current study is based on Advanced Aerodynamic Tools of Large Rotors (AVATAR) European project, and the parameters used for the validation of the computational results are the following ones:

•
Streamwise velocity profiles at different locations (AVATAR Task 3.1 [25])  Turbulence kinetic energy fields at streamwise planes (AVATAR Task 3.2 [26]) Taking into account the previous configuration and in order to characterize the primary vortex generated by the three different VGs, the following parameters have been proposed:

•
Streamwise velocity profiles at different streamwise locations.

•
Primary vortex vertical and lateral positions to define the vortex trajectory.
According to the experimental data available in the AVATAR project and to reproduce the experiments of the VG on a flat plate performed in that project, the numerical simulations have been carried out at a Reynolds number Re = 2600 based on local BL momentum thickness θ = 2.4 mm and with a free stream velocity of U ∞ = 15 m/s.A negligible pressure gradient on the plate is assumed.The angle of incidence of the vane to the oncoming flow is 18 • .The calculation of the local BL momentum thickness θ was made by the application of Equation (1): where u x is the streamwise velocity component, U ∞ the free stream velocity and y the vertical coordinate normal to the wall.In order to validate the numerical model of the VG, the results of the numerical simulations are compared to the experimental data of the AVATAR Project [25,26].The experiments were carried out in the Boundary Layer Wind Tunnel of the Delft University of Technology.The tunnel can attain a maximum speed of 38 m/s in the wide-walled test section of 1.5 × 0.25 m 2 .The large separation of the side walls (1.5 m) minimizes end effects on the flow region of interest along the centerline zone.An adjustable back wall allows adjustment of the pressure gradient, ensuring a truly null pressure gradient when desired.The turbulence level at maximum free stream velocity was determined as 0.5%.Vortex generator dimensions were designed according to the previous research made by m, where the rectangular vortex generators were designed according to the study of Godard et al. [28] and are summarized in Table 1.The placement of the VGs is such that flow impinges on the vanes in a divergent manner, producing counter-rotating, common down-flow embedded vortices.For each test case, 25 pairs of rectangular VGs were mounted in an array configuration, side by side, as shown in Figure 2, covering a spanwise distance of 0.75 m and centered on the tunnel centerline.This was done in order to minimize end-effects from the finite array.The vortex generator height is H = 5 mm and considering a negligible pressure gradient.gradient when desired.The turbulence level at maximum free stream velocity was determined as 0.5%.Vortex generator dimensions were designed according to the previous research made by Baldacchino et al. [27], where the rectangular vortex generators were designed according to the study of Godard et al. [28] and are summarized in Table 1.The placement of the VGs is such that flow impinges on the vanes in a divergent manner, producing counter-rotating, common down-flow embedded vortices.For each test case, 25 pairs of rectangular VGs were mounted in an array configuration, side by side, as shown in Figure 2, covering a spanwise distance of 0.75 m and centered on the tunnel centerline.This was done in order to minimize end-effects from the finite array.The vortex generator height is H = 5 mm and considering a negligible pressure gradient.

Computational Setup
The present study consists of a rectangular VG positioned on a flat plate with an incident angle to the oncoming flow.The height of the VG located at the flat plate in the simulations is going to be different, and the AR will consequently change, whereas the angle of attack will remain constant.The ARs of the three cases are ARH2 = 2, ARH = 2.5 and ARH1 = 3, and the heights are going to be H2 = 4.16, H = 5 and H1 = 6.25 mm, respectively.The computational domain dimensions are represented in Figure 3 normalized by the baseline VG case height H.
Wake symmetry behind the VGs is assumed in the current simulations.Thus, the computational domain includes only one VG inclined to the oncoming flow instead of a pair.In that way, meshing and computational time can be considerably reduced.

Computational Setup
The present study consists of a rectangular VG positioned on a flat plate with an incident angle to the oncoming flow.The height of the VG located at the flat plate in the simulations is going to be different, and the AR will consequently change, whereas the angle of attack will remain constant.The ARs of the three cases are AR H2 = 2, AR H = 2.5 and AR H1 = 3, and the heights are going to be H 2 = 4.16, H = 5 and H 1 = 6.25 mm, respectively.The computational domain dimensions are represented in Figure 3 normalized by the baseline VG case height H.
Wake symmetry behind the VGs is assumed in the current simulations.Thus, the computational domain includes only one VG inclined to the oncoming flow instead of a pair.In that way, meshing and computational time can be considerably reduced.Figure 2 illustrates in green color the spanwise locations of the domain consisting of only one vane.The symmetry assumption used in the present computations can be justified by the previous studies of Sorensen et al. [29] and Velte et al. [12].The computational domain is divided into 28 blocks.The mesh has been refined near the VG and in the corners of the vane where the velocity gradients are large.In regions far away from the VG and the wake, the mesh density is lower.Five blocks are located around the VG and six blocks downstream of the vane to capture the generated primary vortex.The same block-based meshing strategy as in the study of Urkiola et al. [15] has been followed.The total amount of cells is eight million, with a height of the first cell of Δy/H = 3.23 × 10 −6 , normalized according the baseline VG height.Around the vane, the mesh has 1.7 million cells, while the mesh downstream the VG for capturing the wake has approximately 2.4 million cells; see Figure 4.In order to resolve the boundary layer, cell clustering has been used close to the wall, and the wall dimensionless distance of the first layer of cells is less than y+ < 1.Each surface of the domain has been assigned to a type of boundary condition.The vane's four different patches and the bottom of the domain were defined as a wall with a non-slip condition.The roof of the domain and the two lateral surfaces were defined as slip surfaces (symmetrical The computational domain is divided into 28 blocks.The mesh has been refined near the VG and in the corners of the vane where the velocity gradients are large.In regions far away from the VG and the wake, the mesh density is lower.Five blocks are located around the VG and six blocks downstream of the vane to capture the generated primary vortex.The same block-based meshing strategy as in the study of Urkiola et al. [15] has been followed.The total amount of cells is eight million, with a height of the first cell of ∆y/H = 3.23 × 10 −6 , normalized according the baseline VG height.Around the vane, the mesh has 1.7 million cells, while the mesh downstream the VG for capturing the wake has approximately 2.4 million cells; see Figure 4.In order to resolve the boundary layer, cell clustering has been used close to the wall, and the wall dimensionless distance of the first layer of cells is less than y+ < 1.The computational domain is divided into 28 blocks.The mesh has been refined near the VG and in the corners of the vane where the velocity gradients are large.In regions far away from the VG and the wake, the mesh density is lower.Five blocks are located around the VG and six blocks downstream of the vane to capture the generated primary vortex.The same block-based meshing strategy as in the study of Urkiola et al. [15] has been followed.The total amount of cells is eight million, with a height of the first cell of Δy/H = 3.23 × 10 −6 , normalized according the baseline VG height.Around the vane, the mesh has 1.7 million cells, while the mesh downstream the VG for capturing the wake has approximately 2.4 million cells; see Figure 4.In order to resolve the boundary layer, cell clustering has been used close to the wall, and the wall dimensionless distance of the first layer of cells is less than y+ < 1.Each surface of the domain has been assigned to a type of boundary condition.The vane's four different patches and the bottom of the domain were defined as a wall with a non-slip condition.The roof of the domain and the two lateral surfaces were defined as slip surfaces (symmetrical Each surface of the domain has been assigned to a type of boundary condition.The vane's four different patches and the bottom of the domain were defined as a wall with a non-slip condition.The roof of the domain and the two lateral surfaces were defined as slip surfaces (symmetrical hypothesis).The velocity inlet was chosen for the entry of the fluid and pressure outlet for the exit of the fluid downstream the VG.
The quality of the mesh was evaluated by five main indicators to analyze whether it can be classified as a high quality mesh; see Table 2.This set of parameters is a mix of industry standards, solver manuals and academia standards and should consequently not be regarded as the only choice of parameter selection.Equiangular skewness is an indicator of how optimal the cell shape is related to the corner angles.For hexahedral cells, skewness should not exceed 0.85 to obtain a fairly accurate solution.The cell aspect ratio of a cell is typically the ratio between the width and height.For critical flow areas, except those close to the wall, the cell aspect ratio should not exceed an average of 10.  hypothesis).The velocity inlet was chosen for the entry of the fluid and pressure outlet for the exit of the fluid downstream the VG.The quality of the mesh was evaluated by five main indicators to analyze whether it can be classified as a high quality mesh; see Table 2.This set of parameters is a mix of industry standards, solver manuals and academia standards and should consequently not be regarded as the only choice of parameter selection.Equiangular skewness is an indicator of how optimal the cell shape is related to the corner angles.For hexahedral cells, skewness should not exceed 0.85 to obtain a fairly accurate solution.The cell aspect ratio of a cell is typically the ratio between the width and height.For critical flow areas, except those close to the wall, the cell aspect ratio should not exceed an average of 10. Figure 5 represents the wall y+ field distribution on the VG and bottom walls.For the current study, the maximum value of y+ corresponds to 0.366 and its average is 0.044.In this work, the open source code OpenFOAM (Version 2.4.0,The OpenFOAM Foundation Ltd., London, UK) [30] has been used for simulating a rectangular VG on a flat plate with negligible pressure gradient.This CFD code is an object-oriented library written in C++ to solve computational continuum mechanics problems.The solver potential Foam, which solves potential flows, was used to generate starting fields (initialization of fields) in order to speed up the convergence process.This solver is suitable to generate initial conditions for more advanced solvers such as the one used in the present work named simpleFoam.The simpleFoam solver has been applied for steady-state, incompressible and turbulent flows using the RANS (Reynolds averaged Navier-Stokes) equations.Second order discretization schemes were employed in the CFD simulations.The k-omega SST (shear stress transport) turbulence model developed by Menter [31] was used for all the computations.This turbulence model is a combination of two models: Wilcox's k-ω model for the near wall region and the k-ε model for the outer region and in free shear flows.
Simulations were performed on a personal server-clustered parallel machine with Intel ® Core i7-6700 CPU 3.40 GHz × 8 cores and 32 GB RAM on 64-bit Linux.The domain was automatically divided into eight subdomains to solve in parallel and to reduce the simulation time.The computational cost was approximately 12 days per simulation.
A mesh independency study has been done with the boundary layer velocity profile at the plane 10 H downstream of the VG and at z = 0 cross-wise position.Figure 6a shows a comparison between In this work, the open source code OpenFOAM (Version 2.4.0,The OpenFOAM Foundation Ltd., London, UK) [30] has been used for simulating a rectangular VG on a flat plate with negligible pressure gradient.This CFD code is an object-oriented library written in C++ to solve computational continuum mechanics problems.The solver potential Foam, which solves potential flows, was used to generate starting fields (initialization of fields) in order to speed up the convergence process.This solver is suitable to generate initial conditions for more advanced solvers such as the one used in the present work named simpleFoam.The simpleFoam solver has been applied for steady-state, incompressible and turbulent flows using the RANS (Reynolds averaged Navier-Stokes) equations.Second order discretization schemes were employed in the CFD simulations.The k-omega SST (shear stress transport) turbulence model developed by Menter [31] was used for all the computations.This turbulence model is a combination of two models: Wilcox's k-ω model for the near wall region and the k-ε model for the outer region and in free shear flows.
Simulations were performed on a personal server-clustered parallel machine with Intel ® Core i7-6700 CPU 3.40 GHz × 8 cores and 32 GB RAM on 64-bit Linux.The domain was automatically divided into eight subdomains to solve in parallel and to reduce the simulation time.The computational cost was approximately 12 days per simulation.
A mesh independency study has been done with the boundary layer velocity profile at the plane 10 H downstream of the VG and at z = 0 cross-wise position.Figure 6a shows a comparison between the CFD curves of three studied meshes: 2 million (coarse), 4 million (medium) and 8 million (fine) cells for the baseline case with a VG height of H = 5 mm. Figure 6b represents the BL velocity profiles of the fine mesh (blue) versus the experimental one (black) provided by [25].
Appl.Sci.2018, 8, 138 7 of 16 the CFD curves of three studied meshes: 2 million (coarse), 4 million (medium) and 8 million (fine) cells for the baseline case with a VG height of H = 5 mm. Figure 6b represents the BL velocity profiles of the fine mesh (blue) versus the experimental one (black) provided by [25].Verification of sufficient mesh resolution was performed based on the Richardson extrapolation method and applied to the normalized peak vorticity ωx/(U∞/H) at the plane position of 10H behind the vane.Table 3 shows the results of the mesh independency study.A monotonic convergence was achieved.The normalized streamwise velocity curve of the fine mesh and the experimental one presented in Figure 6b fit reasonably well.This higher resolution mesh has been used for the current computations and applied to the three low profile vanes H, H1 and H2.

Results
Numerical simulations of the vortex generating vane on a flat plate were performed.Three dimensions of a single low-profile vane H = 5 mm, H1 = 6.25 mm and H2 = 4.16 mm at an incident angle to the oncoming flow of β = 18° were computed and compared with experimental results.The extraction of the data from the computations was conducted in a similar way to the procedure described in [14], in planes normal to the wall downstream of the VG to capture the development of the vortex.Table 4 shows the height H of the baseline VG and the heights in the other simulated cases H1 and H2.The VG H1 has been designed with an aspect ratio AR = 2, which is the highest one, and the VG H2 with an aspect ratio AR = 3. Verification of sufficient mesh resolution was performed based on the Richardson extrapolation method and applied to the normalized peak vorticity ω x /(U ∞ /H) at the plane position of 10H behind the vane.Table 3 shows the results of the mesh independency study.A monotonic convergence was achieved.The normalized streamwise velocity curve of the fine mesh and the experimental one presented in Figure 6b fit reasonably well.This higher resolution mesh has been used for the current computations and applied to the three low profile vanes H, H 1 and H 2 .

Results
Numerical simulations of the vortex generating vane on a flat plate were performed.Three dimensions of a single low-profile vane H = 5 mm, H 1 = 6.25 mm and H 2 = 4.16 mm at an incident angle to the oncoming flow of β = 18 • were computed and compared with experimental results.The extraction of the data from the computations was conducted in a similar way to the procedure described in [14], in planes normal to the wall downstream of the VG to capture the development of the vortex.Table 4 shows the height H of the baseline VG and the heights in the other simulated cases H 1 and H 2 .The VG H 1 has been designed with an aspect ratio AR = 2, which is the highest one, and the VG H 2 with an aspect ratio AR = 3.First of all, a comparison between experimental and CFD streamwise velocity profiles has been carried out for the baseline case H at the planes 10H, 25H and 50H downstream of the vane; see Figure 7. Additionally, four different cross-wise locations are analyzed, which are at z = 0, z = TE (VG trailing edge), z = D/3 and z = D/2, where D = 30 mm.At the crosswise location of z = 0, the numerical results fit quite well to the experimental ones.As the flow moves away from the center z = 0 to the spanwise direction, some differences can be found between the experimental and the numerical streamwise velocity profiles.First of all, a comparison between experimental and CFD streamwise velocity profiles has been carried out for the baseline case H at the planes 10H, 25H and 50H downstream of the vane; see Figure 7. Additionally, four different cross-wise locations are analyzed, which are at z = 0, z = TE (VG trailing edge), z = D/3 and z = D/2, where D = 30 mm.At the crosswise location of z = 0, the numerical results fit quite well to the experimental ones.As the flow moves away from the center z = 0 to the spanwise direction, some differences can be found between the experimental and the numerical streamwise velocity profiles.According to the information available in the AVATAR project, three variables have been chosen for qualitative comparison between the numerical and the experimental results: the streamwise velocity, the normalized streamwise vorticity and the turbulence kinetic energy.The experimental fields were only available at the downstream plane positions from 10H-50H.According to the information available in the AVATAR project, three variables have been chosen for qualitative comparison between the numerical and the experimental results: the streamwise velocity, the normalized streamwise vorticity and the turbulence kinetic energy.The experimental fields were only available at the downstream plane positions from 10H-50H.
Figure 8 shows the out-of-plane streamwise velocity fields on spanwise planes at three different locations 5H, 10H, 25H and 50H downstream of the vortex generator.The comparison includes the results of the simulations corresponding to the three different VG heights and the experimental data.The height of the vortex generated by the vane increases as the vortex is moving away from the flat plate in the streamwise direction.At the farthest plane position studied 50H, the two vortexes that are visible in the near wake plane positions seem to be unique.The roll-up process largely occurs before the 10H location.The wake evolution is typical of the induced flow field where the velocity deficit is clearly observed in the places nearer the core of the vortex.The vortex size and the velocity at its center increase as it moves away from the flat plate.The velocity in the middle of the vortex increases as it goes away, and the vortex increases, as can be seen by the reduction of the yellow.
The normalized vorticity ωx/(U∞/H) fields at different streamwise planes are illustrated in Figure 9.The vorticity is a pseudovector field, which describes the local spinning motion of a continuum near some point (the tendency of something to rotate), as would be seen by an observer located at that point and traveling along with the flow.As shown in Figure 9, the experimental data and the CFD case of H = 5 mm are very similar, endorsing the accuracy of the current computations.As expected, the vorticity in the core of the vortex decreases as the distance increases, and that decrease is more significant at the planes 25H and 50H.The higher the vane, the larger is the vorticity in the vortex core.The height of the vortex generated by the vane increases as the vortex is moving away from the flat plate in the streamwise direction.At the farthest plane position studied 50H, the two vortexes that are visible in the near wake plane positions seem to be unique.The roll-up process largely occurs before the 10H location.The wake evolution is typical of the induced flow field where the velocity deficit is clearly observed in the places nearer the core of the vortex.The vortex size and the velocity at its center increase as it moves away from the flat plate.The velocity in the middle of the vortex increases as it goes away, and the vortex increases, as can be seen by the reduction of the yellow.
The normalized vorticity ω x /(U ∞ /H) fields at different streamwise planes are illustrated in Figure 9.The vorticity is a pseudovector field, which describes the local spinning motion of a continuum near some point (the tendency of something to rotate), as would be seen by an observer located at that point and traveling along with the flow.As shown in Figure 9, the experimental data and the CFD case of H = 5 mm are very similar, endorsing the accuracy of the current computations.As expected, the vorticity in the core of the vortex decreases as the distance increases, and that decrease is more significant at the planes 25H and 50H.The higher the vane, the larger is the vorticity in the vortex core.The comparison of the turbulence kinetic energy fields is represented in Figure 10.The turbulence kinetic energy (TKE) is the mean kinetic energy per unit mass associated with eddies in turbulent flow and is a measure of the intensity of the turbulence.It is defined as follows: where ux, uy and uz are the velocity field components.As shown in Figure 10, the TKE field comparison between the experimental case and the simulations is very similar at the planes 25H and 50H, whereas at the plane 10H, the experimental case does not have a quick development as the simulation has.On the one hand, for the case of the largest vane height H1, the intensity of turbulence results in a higher y coordinate than the baseline H and then the shortest VG H2.That means that as a consequence of the height of the VG, the turbulent flow moves through at a higher elevation.On the other hand, as the flow moves away from the trailing edge of the actuator, the energy is dissipated to a greater extent in the smallest vane H2.Furthermore, it is clearly visible that for the case of H1 = 6.25 mm at the plane 5H behind the VG, the energy created as a consequence of the vortex is the maximum one in comparison with the other cases H and H2.
Figure 11 shows the normalized streamwise velocity (ux/U∞) profiles for the three different VG heights H1 = 4.16 mm, H = 5 mm and H2 = 6.25 mm at the 10H, 25H and 50H downstream positions.As shown in Figure 11, at the downstream position 25H and spanwise position z = D/2 and at 50H and spanwise positions z = D/2 and z = D/3, the differences of the curves are notable.
In order to analyze the vortex evolution, up to three parameters were previously identified in past studies made by Fernandez-Gamiz et al. [32] and Godard et al. [28] on the streamwise longitudinal vortex embedded in turbulent boundary layers.Therefore, the proposed parameters to characterize the primary vortex generated by a passive rectangular VG are the streamwise peak vorticity ωx,max , the vortex path and the wall shear stress (WSS).The peak vorticity is used to locate the center of the streamwise vortex in each plane position downstream of the VG and, subsequently, to determine the The comparison of the turbulence kinetic energy fields is represented in Figure 10.The turbulence kinetic energy (TKE) is the mean kinetic energy per unit mass associated with eddies in turbulent flow and is a measure of the intensity of the turbulence.It is defined as follows: where u x , u y and u z are the velocity field components.As shown in Figure 10, the TKE field comparison between the experimental case and the simulations is very similar at the planes 25H and 50H, whereas at the plane 10H, the experimental case does not have a quick development as the simulation has.On the one hand, for the case of the largest vane height H 1 , the intensity of turbulence results in a higher y coordinate than the baseline H and then the shortest VG H 2 .That means that as a consequence of the height of the VG, the turbulent flow moves through at a higher elevation.On the other hand, as the flow moves away from the trailing edge of the actuator, the energy is dissipated to a greater extent in the smallest vane H 2 .Furthermore, it is clearly visible that for the case of H 1 = 6.25 mm at the plane 5H behind the VG, the energy created as a consequence of the vortex is the maximum one in comparison with the other cases H and H 2 .
Figure 11 shows the normalized streamwise velocity (u x /U ∞ ) profiles for the three different VG heights H 1 = 4.16 mm, H = 5 mm and H 2 = 6.25 mm at the 10H, 25H and 50H downstream positions.As shown in Figure 11, at the downstream position 25H and spanwise position z = D/2 and at 50H and spanwise positions z = D/2 and z = D/3, the differences of the curves are notable.
In order to analyze the vortex evolution, up to three parameters were previously identified in past studies made by Fernandez-Gamiz et al. [32] and Godard et al. [28] on the streamwise longitudinal vortex embedded in turbulent boundary layers.Therefore, the proposed parameters to characterize the primary vortex generated by a passive rectangular VG are the streamwise peak vorticity ω x,max , the vortex path and the wall shear stress (WSS).The peak vorticity is used to locate the center of the streamwise vortex in each plane position downstream of the VG and, subsequently, to determine the vortex trajectory.In addition, the WSS is a parameter associated with the vortex induced rate of mixing of the outer flow with the boundary layer.The trajectory of the primary vortex generated by a vane-type VG plays a determining role in the mixing performance.vortex trajectory.In addition, the WSS is a parameter associated with the vortex induced rate of mixing of the outer flow with the boundary layer.The trajectory of the primary vortex generated by a vane-type VG plays a determining role in the mixing performance.Figure 12 shows the normalized peak vorticity (ωx,max•H)/U∞ development of the three different heights of the VGs, which in fact represents the decay of the vortex along the streamwise direction.As expected, the vortex decay seems to follow an exponential law in all cases.At the closest positions to the VG, the largest value corresponds to the lowest VG height of 4.16 mm, but that trend is turned around at the distance of 50H in which the highest value is that corresponding to the case with a 6.25-mm vane height.The values of the peak vorticity of the different VG cases almost overlap each other, and consequently, no significant variation was found.
By observing the location of the vortex center as a function of the downstream axis x, one can determine the path of the primary vortex in both lateral (z) and vertical (y) directions.The calculation of the vortex center is made by finding the maximum streamwise vorticity point in the corresponding downstream plane position; see Fernandez-Gamiz et al. [32] and Martinez-Filgueira et al. [21].Figures 13 and 14 show a comparison between the trajectories described by the primary vortex generated by the different VG heights.Streamwise coordinates and both lateral and vertical coordinates are normalized by the baseline VG height H.As illustrated in Figure 13, the lateral path of the three cases tends to follow the direction in which the device is inclined, and the trajectories roughly collapse in the streamwise plane positions studied.The vortex paths in the vertical direction are represented in Figure 14.They collapse in the near wake, from the formation of the vortex up to approximately the downstream plane 18H.As the distance increases far away in the streamwise direction, the higher vortex trajectory corresponds clearly to the higher VG H1 = 6.25.In the far wake positions, a clear gap in the vertical paths is visible between the three VG cases, which does not occur in the lateral paths.Figure 12 shows the normalized peak vorticity (ω x,max •H)/U ∞ development of the three different heights of the VGs, which in fact represents the decay of the vortex along the streamwise direction.As expected, the vortex decay seems to follow an exponential law in all cases.At the closest positions to the VG, the largest value corresponds to the lowest VG height of 4.16 mm, but that trend is turned around at the distance of 50H in which the highest value is that corresponding to the case with a 6.25-mm vane height.The values of the peak vorticity of the different VG cases almost overlap each other, and consequently, no significant variation was found.
By observing the location of the vortex center as a function of the downstream axis x, one can determine the path of the primary vortex in both lateral (z) and vertical (y) directions.The calculation of the vortex center is made by finding the maximum streamwise vorticity point in the corresponding downstream plane position; see Fernandez-Gamiz et al. [32] and Martinez-Filgueira et al. [21].Figures 13 and 14 show a comparison between the trajectories described by the primary vortex generated by the different VG heights.Streamwise coordinates and both lateral and vertical coordinates are normalized by the baseline VG height H.As illustrated in Figure 13, the lateral path of the three cases tends to follow the direction in which the device is inclined, and the trajectories roughly collapse in the streamwise plane positions studied.The vortex paths in the vertical direction are represented in Figure 14.They collapse in the near wake, from the formation of the vortex up to approximately the downstream plane 18H.As the distance increases far away in the streamwise direction, the higher vortex trajectory corresponds clearly to the higher VG H 1 = 6.25.In the far wake positions, a clear gap in the vertical paths is visible between the three VG cases, which does not occur in the lateral paths.The WSS is the result of friction within the fluid and between the fluid and the walls and is related to the fluid viscosity.The magnitude of the WSS is dependent on how fast the velocity increases when moving away from the wall.The wall shear stress, τω, can be determined by Equation (3): where μ is the dynamic viscosity, ux is the axial flow velocity and y is the distance normal to the wall.The units in IS are Pa. Figure 15 represents the comparison of the wall shear stress for the three analyzed VG geometries on a flat plate.Data were extracted from the VG trailing-edge up to the end of the domain and at a spanwise coordinate of z = 0.A dependency between the vane height and the wall shear stress produced in the VG wake is clearly seen in Figure 15.The largest value of the WSS is reached by the highest VG H1 = 6.25 mm.According to Gad-el-Hak [33], transferring momentum towards the near wall region means that the velocity is increased, which leads to an increase in the wall shear stress.In the study of Lin et al. [34] on an airfoil with micro-VGs, it was demonstrated that an increase of the WSS resulted in a lift increasing.The WSS is the result of friction within the fluid and between the fluid and the walls and is related to the fluid viscosity.The magnitude of the WSS is dependent on how fast the velocity increases when moving away from the wall.The wall shear stress, τω, can be determined by Equation (3): where μ is the dynamic viscosity, ux is the axial flow velocity and y is the distance normal to the wall.The units in IS are Pa. Figure 15 represents the comparison of the wall shear stress for the three analyzed VG geometries on a flat plate.Data were extracted from the VG trailing-edge up to the end of the domain and at a spanwise coordinate of z = 0.A dependency between the vane height and the wall shear stress produced in the VG wake is clearly seen in Figure 15.The largest value of the WSS is reached by the highest VG H1 = 6.25 mm.According to Gad-el-Hak [33], transferring momentum towards the near wall region means that the velocity is increased, which leads to an increase in the wall shear stress.In the study of Lin et al. [34] on an airfoil with micro-VGs, it was demonstrated that an increase of the WSS resulted in a lift increasing.The WSS is the result of friction within the fluid and between the fluid and the walls and is related to the fluid viscosity.The magnitude of the WSS is dependent on how fast the velocity increases when moving away from the wall.The wall shear stress, τ ω , can be determined by Equation (3): where µ is the dynamic viscosity, u x is the axial flow velocity and y is the distance normal to the wall.The units in IS are Pa. Figure 15 represents the comparison of the wall shear stress for the three analyzed VG geometries on a flat plate.Data were extracted from the VG trailing-edge up to the end of the domain and at a spanwise coordinate of z = 0.A dependency between the vane height and the wall shear stress produced in the VG wake is clearly seen in Figure 15.The largest value of the WSS is reached by the highest VG H 1 = 6.25 mm.According to Gad-el-Hak [33], transferring momentum towards the near wall region means that the velocity is increased, which leads to an increase in the wall shear stress.In the study of Lin et al. [34] on an airfoil with micro-VGs, it was demonstrated that an increase of the WSS resulted in a lift increasing.

Conclusions
Vortices generated by a passive rectangular vane-type VG on a flat plate have been studied.Numerical simulations at a Reynolds number Re = 2600 based on the local BL momentum thickness θ = 2.4 mm and free stream velocity U∞ = 15 ms −1 have been carried out using the RANS method.
Computational fluid dynamics simulations are able to reproduce the physics of the primary vortex generated by a rectangular passive VG with reasonably good reliability.As numerical methods are the most frequently-used ones applied in optimizing VGs on wind turbine blades, the extensive information presented in the current study can be used as guidance for successful design and implementation of VGs on wind turbine blades and to carry out parametric dependencies of the VGs on different flow conditions.Because of the large range of parameters inherent in the problem (incident angle, interspacing between vanes, device orientation, vane height and length, etc.), an in-depth understanding of the fluid flow is needed to be correctly applied.
Furthermore, a qualitative comparison of the spanwise velocity, normalized vorticity and turbulence kinetic energy fields at different plane positions downstream of the VG has been introduced to see how the vortex is developed.It is verified that when using an appropriate mesh, fields such as axial velocity, vorticity and TKE are determined with reasonable accuracy.
A correlation study and comparison of the vortex parameters such as vortex trajectory, peak vorticity and streamwise BL velocity profiles with three different VG heights has also been carried out.It should be noted that, in the development of the center vortex trajectory, the points overlap each other in the lateral path for the three different VG heights, but not in the vertical path.
It has been demonstrated for the three analyzed geometries that, for the same level of vorticity, the vertical path of the primary vortex and the wall shear stress downstream of the vane are higher as the actuator height increases.The wall shear stress has a direct influence on the quantity of energy transfer.Then, if a delayed stall is desired, the best geometry will be the one with the best wall shear stress increment.Thus, it could be concluded that the highest VG H1 = 6.25 mm is more suitable for separation control applications.Country UPV/EHU through the SAIOTEK (S-PE11UN112) and EHU12/26 research programs, respectively, is gratefully acknowledged.
Author Contributions: Unai Fernandez-Gamiz, Iñigo Errasti and Ruben Gutierrez prepared and ran the numerical part.They also wrote the manuscript.The post-processing was done by Ruben Gutierrez.Ana Boyano and Oscar Barambones provided constructive instructions in the process of preparing the paper.

Conclusions
Vortices generated by a passive rectangular vane-type VG on a flat plate have been studied.Numerical simulations at a Reynolds number Re = 2600 based on the local BL momentum thickness θ = 2.4 mm and free stream velocity U ∞ = 15 ms −1 have been carried out using the RANS method.
Computational fluid dynamics simulations are able to reproduce the physics of the primary vortex generated by a rectangular passive VG with reasonably good reliability.As numerical methods are the most frequently-used ones applied in optimizing VGs on wind turbine blades, the extensive information presented in the current study can be used as guidance for successful design and implementation of VGs on wind turbine blades and to carry out parametric dependencies of the VGs on different flow conditions.Because of the large range of parameters inherent in the problem (incident angle, interspacing between vanes, device orientation, vane height and length, etc.), an in-depth understanding of the fluid flow is needed to be correctly applied.
Furthermore, a qualitative comparison of the spanwise velocity, normalized vorticity and turbulence kinetic energy fields at different plane positions downstream of the VG has been introduced to see how the vortex is developed.It is verified that when using an appropriate mesh, fields such as axial velocity, vorticity and TKE are determined with reasonable accuracy.
A correlation study and comparison of the vortex parameters such as vortex trajectory, peak vorticity and streamwise BL velocity profiles with three different VG heights has also been carried out.It should be noted that, in the development of the center vortex trajectory, the points overlap each other in the lateral path for the three different VG heights, but not in the vertical path.
It has been demonstrated for the three analyzed geometries that, for the same level of vorticity, the vertical path of the primary vortex and the wall shear stress downstream of the vane are higher as the actuator height increases.The wall shear stress has a direct influence on the quantity of energy transfer.Then, if a delayed stall is desired, the best geometry will be the one with the best wall shear stress increment.Thus, it could be concluded that the highest VG H 1 = 6.25 mm is more suitable for separation control applications.
Appl.Sci.2018, 8, 138 3 of 16 two height variations are also investigated, those corresponding to ARH2 = 2 and ARH1 = 3, as sketched in Figure 1.For this purpose, computational simulations have been carried out using the RANS (Reynolds averaged Navier-Stokes) method and at a Reynolds number Re = 2600.The case consists of a single vortex generator on a flat plate with the angle of attack of the vane to the oncoming flow β = 18°.The flow over the flat plate without the VG has been previously simulated to calculate the boundary layer thickness at the VG position.

Figure 1 .
Figure 1.Sketch of the vortex generator (VG) dimensions with respect to the local boundary layer (not to scale).

Figure 1 .
Figure 1.Sketch of the vortex generator (VG) dimensions with respect to the local boundary layer (not to scale).
Figure 2 illustrates in green color the spanwise locations of the domain consisting of only one vane.The symmetry assumption used in the present computations can be justified by the previous studies of Sorensen et al. [29] and Velte et al. [12].

Figure 3 .
Figure 3. Sketch of the computational domain.Dimensions are scaled by the baseline VG height H (not to scale).

Figure 4 .
Figure 4. View of the meshing around the rectangular VG.

Figure 3 .
Figure 3. Sketch of the computational domain.Dimensions are scaled by the baseline VG height H (not to scale).

Figure 3 .
Figure 3. Sketch of the computational domain.Dimensions are scaled by the baseline VG height H (not to scale).

Figure 4 .
Figure 4. View of the meshing around the rectangular VG.

Figure 4 .
Figure 4. View of the meshing around the rectangular VG.

Figure 5
Figure5represents the wall y+ field distribution on the VG and bottom walls.For the current study, the maximum value of y+ corresponds to 0.366 and its average is 0.044.

Figure 5 .
Figure 5. Wall y+ field distribution on the VG and bottom walls.

Figure 5 .
Figure 5. Wall y+ field distribution on the VG and bottom walls.

Figure 6 .
Figure 6.Comparison of the normalized streamwise velocity profile at the 10H downstream plane position and spanwise position of z = 0 corresponding to the baseline case of H = 5 mm.(a) CFD curves for three different meshes: coarse (green), medium (red) and fine (blue); (b) Boundary layer (BL) velocity profiles of the fine mesh (blue) versus the experimental one (black).

Figure 6 .
Figure 6.Comparison of the normalized streamwise velocity profile at the 10H downstream plane position and spanwise position of z = 0 corresponding to the baseline case of H = 5 mm.(a) CFD curves for three different meshes: coarse (green), medium (red) and fine (blue); (b) Boundary layer (BL) velocity profiles of the fine mesh (blue) versus the experimental one (black).

Figure 7 .
Figure 7. Normalized streamwise velocity profiles () at different streamwise locations 10H (top), 25H (middle) and 50H (bottom) behind the VG.Solid black lines (-) correspond to the experimental data and the blue ones (-) to numerical profiles.

Figure 7 .
Figure 7. Normalized streamwise velocity profiles () at different streamwise locations 10H (top), 25H (middle) and 50H (bottom) behind the VG.Solid black lines (-) correspond to the experimental data and the blue ones (-) to numerical profiles.

Figure 8
Figure 8 shows the out-of-plane streamwise velocity fields on spanwise planes at three different locations 5H, 10H, 25H and 50H downstream of the vortex generator.The comparison includes the results of the simulations corresponding to the three different VG heights and the experimental data.

Figure 8 .
Figure 8.Comparison of the out-of-plane streamwise velocity fields ux between experimental and numerical simulations of the three cases H, H1 and H2.Streamwise plane positions from top to bottom: 5H, 10H, 25H, 50H.

Figure 8 .
Figure 8.Comparison of the out-of-plane streamwise velocity fields u x between experimental and numerical simulations of the three cases H, H 1 and H 2 .Streamwise plane positions from top to bottom: 5H, 10H, 25H, 50H.

Figure 12 .
Figure 12.The normalized peak vorticity development for the different VG heights.

Figure 12 .
Figure 12.The normalized peak vorticity development for the different VG heights.

Figure 12 .
Figure 12.The normalized peak vorticity development for the different VG heights.

Figure 13 .
Figure 13.The development of the vortex center in the horizontal direction.

Figure 14 .
Figure 14.The development of the vortex center in the vertical direction.

Figure 13 .
Figure 13.The development of the vortex center in the horizontal direction.

Figure 13 .
Figure 13.The development of the vortex center in the horizontal direction.

Figure 14 .
Figure 14.The development of the vortex center in the vertical direction.

Figure 14 .
Figure 14.The development of the vortex center in the vertical direction.

Figure 15 .
Figure 15.Wall shear stress comparison in the streamwise direction from the VG trailing-edge to 50H (spanwise distance z = 0) of the three different VG heights and without any vane on the wall.

Figure 15 .
Figure 15.Wall shear stress comparison in the streamwise direction from the VG trailing-edge to 50H (spanwise distance z = 0) of the three different VG heights and without any vane on the wall.

Table 1 .
Vortex generator geometry in the Advanced Aerodynamic Tools of Large Rotors (AVATAR) project.

Table 1 .
Vortex generator geometry in the Advanced Aerodynamic Tools of Large Rotors (AVATAR) project.

Table 3 .
Results of the mesh independency study 1 .
1RE represents the extrapolated solution, R the ratio of errors and p the order of accuracy.

Table 3 .
Results of the mesh independency study 1 .
1 RE represents the extrapolated solution, R the ratio of errors and p the order of accuracy.