Numerical Study of the Effect of the Reynolds Number and the Turbulence Intensity on the Performance of the NACA 0018 Airfoil at the Low Reynolds Number Regime

Abstract

In recent years, there has been an increased interest in the old NACA four-digit series when designing wind turbines or small aircraft. One of the airfoils frequently used for this purpose is the NACA 0018 profile. However, since 1933, for over 70 years, almost no new experimental studies of this profile have been carried out to investigate its performance in the regime of small and medium Reynolds numbers as well as for various turbulence parameters. This paper discusses the effect of the Reynolds number and the turbulence intensity on the lift and drag coefficients of the NACA 0018 airfoil under the low Reynolds number regime. The research was carried out for the range of Reynolds numbers from 50,000 to 200,000 and for the range of turbulence intensity on the airfoil from 0.01% to 0.5%. Moreover, the tests were carried out for the range of angles of attack from 0 to 10 degrees. The uncalibrated $\gamma -R{e}_{\theta }$ transition turbulence model was used for the analysis. Our research has shown that airfoil performance is largely dependent on the Reynolds number and less on the turbulence intensity. For this range of Reynolds numbers, the characteristic of the lift coefficient is not linear and cannot be analyzed using a single aerodynamic derivative as for large Reynolds numbers. The largest differences in both aerodynamic coefficients are observed for the Reynolds number of 50,000.

1. Introduction

Nowadays, computational fluid dynamics (CFD) methods are still not fast enough and sufficiently precise in terms of large angles of attack to completely replace simplified engineering methods based on the principle of conservation of momentum used to simulate the movement of micro-air vehicles, high-altitude unmanned air vehicles, or rotors of wind turbines [1,2,3,4,5]. The accuracy of the results of the aerodynamic characteristics of both aircraft and rotors obtained by these engineering methods depends largely on the airfoil characteristics [6]. These airfoil characteristics, in the form of the coefficients of the aerodynamic forces as a function of the angle of attack, constitute the input data for the engineering methods [7,8,9,10,11]. In recent years, there has been an increased interest in small aircraft and small wind turbines [12,13,14,15,16,17]. From the point of view of simplified approaches, a significant problem is the working conditions in which these devices must operate. One of the key quantities determining these working conditions is the Reynolds number and turbulence. The Reynolds number experienced by the airfoils of small wind turbines or UAVs can range from tens of thousands to several hundred thousand [18,19,20,21,22]. Therefore, these are not conditions encountered in traditional aviation. However, in the case of turbulence, the problem is the different turbulence intensity values found in wind tunnels and real operating conditions. Generally, a value of less than 1% is considered to be of low turbulence intensity, while a value of more than 10% is considered to be of high turbulence intensity. In the case of modern low turbulence wind tunnels, the turbulence intensity of undisturbed flow can be even less than 0.05% [23,24]. The turbulence intensity measured in the atmosphere is usually from a few to several dozen percent [25,26]. The aerodynamics of a vertical axis wind turbine rotor (Figure 1) is even more complicated as its blades constantly operate under varying flow conditions. As the rotor rotates, the angle of attack of the airfoil changes, and thus the local Reynolds number also changes constantly [27,28,29,30,31]. During the operation of the rotor of the Darrieus wind turbine, the rotor blade works with variable turbulence intensities [32]. This is due to the slowdown and disturbance of the airflow by a blade moving in the upwind part of the rotor. In the downwind part of the rotor, the blade already works in a sometimes highly disturbed flow.

In the previous work, Rogowski et al. [6] studied the aerodynamic characteristics of the NACA0018 airfoil operating at the low Reynolds number regime ($Re=160\times {10}^3$) and at very low inlet turbulence intensity corresponding to a very clean flow. These authors used the uncalibrated Transition SST approach implemented in the ANSYS Fluent CFD solver. This paper discusses the performance of this profile considering Reynolds numbers from 50,000 to 200,000 and turbulence intensities at the airfoil ranging from 0.01% to 0.5%. The numerical results of the aerodynamic force coefficients of the symmetrical NACA 0018 airfoil shown in this paper are primarily intended to serve as input data for analytical approaches dedicated to Darrieus wind turbine rotors [7,33].

The result of the work of the American agency NACA (National Advisory Committee for Aeronautics) was the preparation of Technical Report No. 460 “The Characteristics of 78 Related Airfoil Sections from Tests in the Variable-Density Wind Tunnel” [34]. It presents a new family of four-digit airfoils intended to operate in the high Reynolds number regime. Researchers found that the airfoil performance is most influenced by three parameters: the maximum camber, the location of the maximum camber, and the maximum thickness of the airfoil. The NACA abcd profile designation contains information on these three parameters: the first digit a denotes the maximum camber as a percentage of the chord, the second digit b indicates the distance from the leading edge to the maximum camber in tenths of the chord, whereas the last two digits cd describe the maximum thickness of the profile as a percent of the chord [35]. For example, in the case of the NACA 0018 profile, the maximum camber is 0%, and the maximum thickness is 18%. The airfoil’s maximum thickness is located at 30% of the chord length from the leading edge. The airfoil shape is described using a mathematical formula available in many scientific publications, e.g., [36,37]. Even though the NACA 0018 airfoil was designed relatively long ago, it is still prevalent in VAWT design because of the balance between structural integrity and aerodynamic efficiency [24,38].

The famous Report No. 460 from 1933 published the NACA 0018 profile characteristics for a Reynolds number of 3,200,000 [34]. However, these results by the American agency NACA are not helpful as inputs for engineering models developed to find the aerodynamic characteristics of Darrieus wind turbines because the blade Reynolds numbers of even the largest wind turbines of this type are much lower [39,40]. More valuable characteristics can be found in Report No. 568 published in 1937, which presents the NACA 0018 airfoil performance measured for the Reynolds number regime from 40,000 to 3,000,000 [41]. The characteristics of the lift coefficients presented in this report show some, though slight, dependence on the Reynolds number. A much more significant effect of the Reynolds number was observed for the drag coefficient. In the case of the largest Reynolds numbers analyzed by these authors ($2.97\times {10}^6$ and $2.36\times {10}^6$), an almost independence of the lift coefficients from the Reynolds number in the linear part of the $C_L\left(\alpha \right)$ curve is observed, and a slight influence of the Reynolds number on the drag coefficient [41]. A similar trend of the $C_L$ characteristics for the NACA 0009 airfoil for Reynolds numbers in the range of $3\times {10}^6$ to $6\times {10}^6$ can be found in the well-known report “Summary of Airfoil Data” from 1945 [42]. Even though the Jacobs and Sherman results are more useful for VAWTs because of the low Reynolds numbers [41]. Unfortunately, the range of angles of attack considered in this experiment is insufficient to calculate the performance of a Darrieus rotor operating at low tip speed ratios. For this reason, as part of the program on the Darrieus vertical axis wind turbine (VAWT) at Sandia National Laboratories, a series of measurements in the wind tunnel for four airfoil sections: NACA-0009, -0012, -0012H, and -0015 were carried out at low Reynolds numbers. These obtained experimental data of the airfoil sections were then extrapolated for the Reynolds number range from ${10}^4$ to ${10}^7$, and for the NACA-0018, -0021 and -0025 airfoil sections [7,43]. For many years, researchers and engineers have used this database to analyze and design vertical axis wind turbines [44,45,46,47,48]. However, since these NACA 0018 airfoil characteristics obtained by Sheldahl and Klimas were extrapolated (not measured directly), until the beginning of this century, only measurement data published by Jacobs and Sherman in 1937 were available. Of course, there are sometimes individual studies, e.g., [27], of this airfoil in the literature, but they are not carried out for a wide range of flow parameters. Only in the last two decades, as a result of increased interest in vertical-axis wind turbines, several experimental studies of the NACA 0018 profile were carried out [24,49,50,51,52,53,54,55,56]. The relatively late appearance of a new measurement series for the NACA 0018 airfoil for a long time retained the possibility of validating the aerodynamic characteristics calculated using modern analytical tools offering the possibility of capturing transition phenomena. In 2019, Melani et al. [57] published an important paper summarizing the collected results of experimental and numerical studies of this airfoil for low and medium Reynolds numbers. A juxtaposition of the results of various researchers shows significant differences in the lift characteristics for Reynolds numbers below 300,000. This is due to, among other things, the sensitivity to the flow conditions in the wind tunnel. Timmer concluded that the characteristics published in 1937 by Jacobs and Sherman might be affected by errors resulting from the high level of turbulence in the wind tunnel, which may significantly affect the transition characteristics [24,41]. Research by Melani et al. [57] also shows that the $\gamma -R{e}_{\theta }$ turbulence model provides sufficiently accurate aerodynamic characteristics for the NACA 0018 airfoil. As mentioned above, in the case of medium and high Reynolds numbers, its influence on the airfoil performance is not so significant, which is also confirmed by the observations of Melani et al. [57], and sufficient accuracy of the results may also be ensured by, e.g., the popular two-equation k-ω SST turbulence model [20]. In addition to the dependence of the airfoil characteristics on the angle of attack and the Reynolds number, a dependence on the freestream turbulence intensity is also necessary to find the reliable aerodynamic performance of the VAWTs [58]. This is evidenced by some recent precise experimental and numerical studies [53,59]. However, the literature still lacks the NACA 0018 airfoil characteristics calculated with the uncalibrated $\gamma -R{e}_{\theta }$ model for a wide range of different flow parameters to serve as a benchmark for further investigations.

As discussed in the first paragraph of this section, the accuracy of engineering approaches for vertical axis wind turbines depends, inter alia, on the accuracy of the input data of the airfoil characteristics. In aerospace engineering and wind energy, numerical techniques are very often used in parametric studies of airfoils. It is believed that solving the Navier–Stokes (NS) equations will provide information on all turbulence scales in the flow, making it possible to study all the aerodynamic phenomena in a given flow. The Direct Numerical Simulation (DNS) method used, among other things, by Marxen et al. [60] to study the flow over a flat plate allowed a more detailed analysis of the separation-induced transition phenomenon. Another slightly less numerically expensive technique is an approach called large eddy simulation (LES). In the 1960s, this method was proposed to simulate atmospheric flows. With the increase in computing power of supercomputers, this approach has become very attractive and promising for simulating turbulent flows [61,62,63]. From an engineering point of view, both of these methods, DNS and LES, are still numerically too expensive to simulate the wall-bounded flows. In aerospace engineering, Reynolds-averaged Navier–Stokes (RANS) techniques together with one- or two-equation turbulence models are still common [64]. However, numerous studies have shown that classical turbulence models fail for small Reynolds numbers [28,57]. In order to avoid these shortcomings, in recent years, RANS techniques have been developed to capture the separation-induced transition. The $\kappa -{\kappa }_L-\omega$ [65,66] and $\gamma -R{e}_{\theta }$ [67] models are some of the better-known turbulence models available in many modern CFD codes, e.g., [23]. These turbulence models have been successfully validated against standard test cases [65,68]. Suluksna and Juntasaro showed that the predictive reliability of the $\gamma -R{e}_{\theta }$ model was better compared to other correlation-based models [69]. Choudhry and Arjomandi compared these transition models for the NACA 0021 airfoil operating in the small Reynolds number regime [70]. In the conclusions, these authors indicated a slightly higher accuracy of the $\kappa -{\kappa }_L-\omega$ turbulence model in comparison with the $\gamma -R{e}_{\theta }$ model. In recent years, a comparison of five turbulence models has also been made by Aftab et al. [71] for the flow around the NACA 4415 airfoil at a Reynolds number of 120,000. They compared the following approaches: the one-equation Spalart–Allmaras, the two-equation SST k-ω, the SST k-ω with the intermittency equation, the $\kappa -{\kappa }_L-\omega$, and the $\gamma -R{e}_{\theta }$. In the conclusions, the authors indicated reliable results of the aerodynamic force coefficients of both transition models. However, they reported problems with the computation time for the $\kappa -{\kappa }_L-\omega$ model and problems with the accurate capture of the reattachment of the flow for a moderate angle of attack equal to 6 degrees [71]. The CFD-based transition modeling approaches mentioned in this section are not the only ones available on the market. It is also worth mentioning that there are relatively fast analytical tools, such as XFOIL developed by Drela or RFOIL [24], modified versions of the well-known XFOIL code, which are relatively good at determining aerodynamic profile characteristics, taking into account the phenomena related to laminar-turbulent transition. However, these approaches are challenging to implement in modern CFD codes [67]. It should be added that they can be an additional comparative tool for profile performance studies [6,72].

Given the differing opinions on the robustness of the $\kappa -{\kappa }_L-\omega$ approach, the authors of this paper aim to thoroughly investigate the $\gamma -R{e}_{\theta }$ turbulence model available in the ANSYS Fluent CFD code called Transition SST without correcting the default constants of the model [23]. The intention of the authors is, among other things, to provide guidelines for users of this commonly used CFD code as well as for users of other CFD codes with an implementation of the $\gamma -R{e}_{\theta }$ turbulence model and, in addition, to obtain a database for further comparisons.

2. Numerical Model and Its Validation

Rogowski et al. [6] analyzed the well-known NACA 0018 airfoil for only a Reynolds number of 160,000. This was a direct result of the operating conditions of the vertical-axis wind turbine analyzed at TU Delft [73]. Moreover, in the previous work, the two-dimensional airfoil section with chord length c of 1 m was analyzed in a very clean flow—the turbulence intensity was assumed to be 0.05%, which corresponded to the flow conditions in the experiment of Timmer at TU Delft [24]. This paper aims to estimate the aerodynamic performance of the NACA 0018 airfoil at different values of turbulence intensity in undisturbed flow and at different Reynolds numbers.

2.1. Numerical Domain and CFD Solver Settings

As shown in Figure 2, simulations were performed using a C-mesh surrounding the NACA 0018 profile and extending 7.5 c from the semicircular edge of the domain and 15 c from the vertical edge of the domain. The size of the computational domain is vital for the results of the aerodynamic forces. Therefore, this effect will be discussed in Section 2.2 of this paper. The “pressure outlet” boundary condition was assumed on the vertical edge of the domain and the “velocity inlet” boundary condition on the remaining outer edges of the computational domain. The airfoil surfaces are treated as no-slip walls—zero velocity was ensured on the airfoil surfaces.

The two-dimensional, incompressible governing equations are considered in all analyzes presented in this paper. The unsteady Reynolds-averaged Navier–Stokes (URANS) approach with the SIMPLE pressure velocity coupling scheme and second-order discretization for all variables were employed. Moreover, in all simulations presented herein, the Green–Gauss node-based method was employed [74]. The flow over the NACA 0018 profile was investigated using the commercial computational fluid dynamics (CFD) code ANSYS Fluent 19.1 [23].

In this work, the transient flow around a symmetrical airfoil is discussed. The effect of such simulations are, among others, oscillations of aerodynamic forces that result from pressure changes around the airfoil. Observing these transient changes in the airfoil’s aerodynamic characteristics is possible if the time step size is sufficiently small, adapted to the observed physical phenomena occurring during the flow around the airfoil. As shown in the previous research [6], setting even a very small time step size does not give significant changes in the characteristics of the aerodynamic force coefficients. Of course, the oscillation size depends on the angle of attack. As the angle grows, the oscillations become larger and larger. However, even in the case of the largest analyzed angle of attack, equal to 10 degrees, the deviation of the instantaneous values of the aerodynamic force components from the average was not too large. This means that the flow is almost without separation for angles of attack up to 10 degrees. It also means that the instantaneous pressure coefficient distributions resemble those calculated with the Reynolds-averaged Navier–Stokes (RANS) technique. Based on the previous research, in all of the results presented in this paper, the time step size was assumed to be 0.0001 s. For the transient solution using the implicit formulation, ANSYS Fluent CFD code requires multiple iterations to be specified at each time step. In this work, the convergence of every time step solution was usually achieved at approximately 10–15 iterations, with a declaration of 20 maximum iterations per time step. The level of residuals of 10⁻⁶ was used for all computed variables.

An important problem in time-dependent analysis is the effect of initial conditions. In the initial phase of the simulation, the oscillations of the aerodynamic forces are significant, and the values of these forces differ significantly from the final ones. In this work, uniform velocity in all finite volumes of the numerical mesh was assumed as initial conditions. The initial conditions for velocity and turbulence were taken from the boundary conditions at the inlet. Simulations of ten seconds of air movement showed that this effect disappears after approximately 5–6 s. The aerodynamic force coefficients were averaged every second. The aerodynamic characteristics shown in the graphs in the result section of this paper (Section 3) correspond to the last tenth second of motion.

2.2. Computing Domain Size Effect

As mentioned above, the size of the computing domain is an essential aspect of CFD calculations. Some authors pointed out that the influence of the domain size on the accuracy of the calculated airfoil characteristics depends, inter alia, on the type of mesh (e.g., C-mesh, O-mesh, etc.), its density, and physical parameters of the flow. A review of the literature shows that there is no consensus on the size of the computational domain concerning the length of the chord. This applies to both airfoils operating at low and high Reynolds numbers [75,76,77,78]. Some researchers employ small computing domain sizes. In this case, the distance from the obstacle to the domain’s boundary is from several to dozen airfoil chords. Conventionally, the domain area in the wake region is longer than the distance from the inlet to the obstacle. Some researchers create computing domains of huge sizes. Sometimes the distance from the leading edge of the airfoil to the boundary of the computational domain is even several hundred times larger than the airfoil chord length.

Paper [6] discusses the effect of three domain widths, 3.75 c, 7.5 c, and 15 c, on the characteristics of the lift and drag coefficients. Those authors concluded that the domain width does not significantly affect the airfoil aerodynamic characteristics. However, the drag coefficient is more sensitive to changing the domain width than the lift coefficient. A numerical domain with a width of 15 c gives acceptable results for the aerodynamic components. At the same time, the number of mesh nodes is half as compared to the mesh for a domain twice as wide. Therefore, a nominal 15 c-wide domain was used in these studies (the radius of the semicircular outer edge of the domain is 7.5 c).

The length of the computational domain in the wake area behind the obstacle is also an important aspect of CFD calculations. This effect was not discussed in previous studies by Rogowski et al. [6]. Therefore, the authors of this paper decided to verify the numerical model for several different lengths of the computational domain. For this purpose, four computational domains of different lengths were developed. The variable in these analyzes is the length L, measured from the airfoil’s trailing edge to the domain’s outlet boundary, as shown in Figure 3a. The distance from the leading edge of the airfoil to the inlet was the same for all analyzed test cases. Since the range of angles of attack from 0 to 10 degrees is discussed in this paper, the numerical model was verified for two angles of attack: 4 and 10 degrees. Thus, the analysis of the effect of length L on the aerodynamic characteristics of the NACA 0018 airfoil was carried out based on eight test cases. The change in the length of the computational domain results in a change in the number of finite volume elements of the model (the size of the computational grid). For the lowest domain considered in this paper, the number of mesh elements was 797,000, while for the largest domain it was 871,000. Twenty seconds of air movement was simulated for each test case. In addition, a Reynolds number of 150,000 and a turbulence intensity of 0.25% were used to verify the numerical model. Taking into account both the range of Reynolds numbers as well as the range of turbulence intensity analyzed in this paper, the values of Re and TI used for this verification are more or less average values.

The results of the numerical calculations of the lift and drag coefficients are presented in Figure 3b,c, respectively. These charts were made in such a way that it was possible to compare the coefficients both qualitatively and quantitatively. Therefore, the results for the angle of attack equal to 4 degrees are presented on the ordinate axis, whereas the results for the angle of 10 degrees are plotted on the auxiliary y-axis. The results compiled in this way show that the numerical values of the coefficients differ very little for a given angle of attack. For an angle of attack of 4 degrees and a length L of 14 c, the drag coefficient is 0.28% less than the drag coefficient for the case L = 10 c. Figure 3c shows that, in contrast to the drag coefficient, as the length of the computational domain in the wake region increases, the lift coefficient increases. However, this increase is also slight. For an angle of attack of 4 degrees and a length L of 14 c, it is only 0.26% compared to the case of L = 10 c. For the angle of attack of 4 degrees, the increase in the $C_L/C_D$ ratio for the case of L = 14 c is also only 0.54% in relation to the case of L = 10 c, a further increase in the parameter L from 14 c to 20 c causes an increase in the $C_L/C_D$ ratio by only 0.11%. The conducted analyzes confirm that for the L of 14 c, the characteristics of the aerodynamic force coefficients can be treated as largely independent of the domain length in the wake area.

This study also investigated the effect of the computational domain length on the wall Y + parameter and the decrease in turbulence intensity in the area from the inlet to the obstacle. However, the conducted analysis showed that the influence of the computational domain length on both parameters was so negligible that the discussion was limited only to this comment. Both of these parameters are affected, as in the case of aerodynamic forces, only by the angle of attack.

2.3. Turbulence Intensity Effect and Mesh Sensitivity Studies

One of the key achievements of this paper is the analysis of the aerodynamic performance of the symmetrical NACA 0018 airfoil in relation to the turbulence intensity. As a result of numerical dissipation, the turbulence intensity just before the airfoil is sometimes completely different from the turbulence intensity assumed for the boundary conditions. This analysis aims to find a correlation between the turbulence intensity at an inlet and just before the airfoil (at the point located close to the airfoil nose). The different types of boundary conditions available in CFD codes require different flow parameters to be specified. In this work, the velocity inlet boundary condition is employed as shown in Figure 2a. In the current investigations, two turbulence quantities were specified at the inlet: the turbulence intensity, $T{I}_0$, and the turbulence length scale, l. The turbulence intensity is determined as the ratio of the root-mean-square of the velocity fluctuations to the average flow velocity. The turbulence length scale is the quantity that is appropriate to the size of large eddies which contain energy in turbulent flow [23]. Many reports from experimental measurements in wind tunnels only provide information about turbulence intensity at the location of the analyzed obstacle. It is very often the only turbulence parameter measured during experiments. Following the guidelines available in the ANSYS Fluent documentation [23], in the case of wall-bounded flows, it is recommended to calculate the turbulence scale from the following formula: $l=0.4 {\delta }_{99}$, where ${\delta }_{99}$ is the boundary-layer thickness. In the current research, the turbulence length scale, l, was estimated on the airfoil for the coordinate x/c = 0.3, and it is equal to 0.003845 m. This paper does not analyze the sensitivity of the airfoil performance due to this parameter.

A numerical experiment was carried out to investigate the rate of decay of turbulence intensity and determine the correlation between the turbulence intensity at the domain inlet and just before the airfoil. The analysis was performed for the following flow parameters: the Reynolds number was assumed to be 150,000, the undisturbed flow velocity was 2.1911 m/s, and the angle of attack was 6 degrees. In this study, eight values of the turbulence inlet intensity were used: 0.1%, 0.5%, 1%, 2%, 4%, 6%, 8%, and 10%. Unsteady calculations were performed for ten seconds of air movement, assuming a time step length of 0.0001 s. Figure 4b,c illustrates the dependence of the turbulence intensity on the x coordinate, starting from the boundary of the computational domain up to the checkpoint located just at the airfoil nose. Both points representing the start and end of the control segment are shown in Figure 4a. The decay rate of this quantity is strongly dependent on the turbulence intensity that was set for the boundary condition $T{I}_0$. The higher the $T{I}_0$, the higher the decay rate. As shown in Figure 4b, in the case of the lowest turbulence intensity analyzed in this work, equal to 0.1%, the decay rate is represented by an almost constant function slightly dependent on the x coordinate. In the case of turbulence intensity higher than 0.5%, the characteristic inflection of the TI curve is visible. This inflection occurs approximately 0.5 c from the first control point located at the domain entrance. Before this inflection, the turbulence intensity decay is much higher than after it. Near the nose of the airfoil, the turbulence intensity is very low for each test case. Therefore, to compare the differences in the turbulence intensity values in the vicinity of the airfoil, an additional chart was prepared in which the x coordinate was limited to the last meter before the airfoil (Figure 4c). In this figure, it can be seen that the TI just before the airfoil varies from 0.0484% for $T{I}_0=0.1\%$ to $T{I}_0=10\%$. In the first case, it gives a decrease of 51.6%, in the second case by 95.3%. This data can be used to interpolate the $T{I}_0$ value to obtain the expected TI value on the profile. This data shown in Figure 4b,c can be used for interpolation of $T{I}_0$ to obtain a specific TI on the airfoil.

Moreover, a mesh density influence on turbulence decay was examined. Four additional meshes were prepared, from extra coarse to extra fine. Cases with typical settings, Reynolds number, 150,000, and inlet $T{I}_0$ equals 0.25%, were prepared. Studies were performed for two angles of attack: 4 degrees and 10 degrees. All results presented in Table 1 show that differences from the “Medium” case are tiny. Therefore, it can be claimed that mesh density influence on turbulence decay is negligible.

In addition to the rate of turbulence intensity decay with distance from the nose of the airfoil, this paper also discusses the effect of turbulence intensity on the aerodynamic characteristics of the airfoil. These characteristics are discussed in detail in the Results section of this paper. However, Figure 4d–g shows the time dependence of the components of the aerodynamic forces. These graphs show that the coefficients become time independent after 8–9 s. However, the lift coefficient becomes time independent faster than the drag coefficient.

2.4. Validation of the Numerical Model

Figure 5 shows a comparison of the aerodynamic force characteristics of the NACA 0018 airfoil obtained using the $\gamma -R{e}_{\theta }$ approach with the aerodynamic airfoil characteristics measured in the Delft University Low-speed Wind Tunnel by Timmer [24]. The lowest Reynolds number published in the paper [24], equal to $0.15\times {10}^6$, was used to compare the numerical results. With this Reynolds number, the turbulence intensity in the test section was 0.02%. This turbulence intensity value is very close to the lowest one considered in these investigations. The graphs shown in Figure 5 have the same scales as those in Figures 12 and 14 in the Results section of this paper. As shown in Figure 5b, the drag coefficient results agree with the experiment both qualitatively and quantitatively. The absolute error of the minimum drag coefficient compared to the experiment is only 0.00282. In the case of the lift coefficient (Figure 5a), both series of results are sufficiently consistent in the range of angles of attack up to 6 degrees. Slightly more significant differences in the coefficients are visible above this angle of attack. This effect is shown in the earlier work of Rogowski et al. [6] as well as in the studies by Melani et al. [57]. One of the possible causes of differences is the 3D effects not considered in these studies. The second reason is that the turbulence model is not calibrated through the model constants that control the characteristics of the separation bubble [70].

3. Results

3.1. Reynolds Number Effect on Lift and Drag Coefficients

At a low Reynolds number regime, relatively small differences in this quantity can significantly change lift and drag airfoil characteristics. This section details the effect of the Reynolds number on the shape of the $C_L\left(\alpha \right)$ and $C_D\left(\alpha \right)$ curves. It was decided to discuss both components separately in this section.

3.1.1. Lift Force Coefficient

All simulations carried out in this paper, regardless of the Reynolds number and the turbulent intensity showed a zero-lift coefficient at a zero angle of attack. The characteristics of the lift coefficient as a function of the angle of attack at the low Reynolds number regime are slightly different from the characteristics at much higher velocities [41]. As the Reynolds number increases, the influence of the laminar separation bubble on the lift coefficient $C_L$ characteristics decrease, as can be seen in Figure 6.

Analyzing these lift force coefficient characteristics (Figure 6), the characteristic “bend” of the $C_L$ curve can be observed. The lift force characteristics are not linear in the whole range of the angles of attack up to the critical angle. In other words, the derivative $d{C}_L/d\alpha$ is not constant over this range of angles of attack. On both sides of this bend, the $C_L$ characteristics can be analyzed by considering two derivatives $d{C}_L/d\alpha$ that can be regarded to a large extent as constant. Both Gerakopulos et al. [79] and Rogowski et al. [6] used such derivatives to define the two regions. In both of these regions, the increase in the lift coefficient with the angle of attack is largely linear. The region in the range of the angles of attack from zero to the location of the bend is called the “first region”, whereas the region above this bend is called the “second region”. The approach for determining the aerodynamic derivative of $d{C}_L/d\alpha$ used in this work was taken from publication [79]. Figure 7 shows an example of determining these two derivatives of the $C_L$ curve for the Reynolds number of $150\times {10}^3$ and for the turbulence intensity, TI, of 0.25%.

As mentioned above, Figure 6 and Figure 7 only show the lift coefficient results obtained for one value of turbulence intensity of 0.25%. This is due to two reasons. First, this is the average turbulence intensity value for the turbulence intensity range considered in this work. Secondly, the $C_L$ results for the remaining turbulence intensities considered in this paper are very similar. Thirdly, we found the experimental data of the aerodynamic derivatives developed for similar flow conditions. According to Gerakopulos et al. [79], the turbulence intensity measured during the experiment was less than 0.3%. Therefore, there was a possibility of further validation of our results. Figure 8 shows a comparison of the experimental and numerical results of the $d{C}_L/d\alpha$ derivatives separately for two regions. According to the authors of this paper, the accuracy of the numerical series is satisfactory. For the Reynolds number of 100,000, the most significant differences in the results of the $d{C}_L/d\alpha$ derivatives are visible. However, they may result from both the dispersion of the values of the lift coefficients used to calculate the derivatives based on the linear approximation and the range of angles of attack used to calculate individual derivatives.

From the observation of the data shown in Figure 8, the derivatives in the first region tend to decrease as the Reynolds number increases, while the derivatives in the second region are not as dependent on the Reynolds number. As seen “with the naked eye” in this figure, an increase in the Reynolds number causes an increase in the first region. However, a more detailed analysis of this effect requires the location of the transition point between the first and second regions to be specified. As mentioned above, the $C_L\left(\alpha \right)$ is a continuous function around the bend. This publication proposes another way to find the transition location. For this paper, the transition angle of attack was defined and marked with the symbol ${\alpha }_t$. This angle was found analytically by solving a simple system of equations composed of trend line equations approximating both regions. In other words, it is the angle at the intersection of the trend lines for the first and second regions (Figure 9a). Figure 9b quantitatively documents the nearly linear growth of the first region as the Reynolds number increases.

3.1.2. Drag Force Coefficient

Figure 10 shows the dependence of the drag coefficient on the angle of attack and the Reynolds number. These results are also discussed for one turbulence intensity of 0.25%. The largest drag coefficients are calculated for the lowest Reynolds number taken into account in this study, of 50,000. At this Reynolds number, a faster increase in drag is also observed starting from an angle of attack of 6 deg. In the case of Reynolds numbers larger than 50,000, the increase in drag as a function of the angle of attack is much smoother.

In this work, the qualitative dependence of the drag coefficient as a function of the Reynolds number for individual angles of attack was also examined. Figure 11 shows the dependence of the drag coefficient on the Reynolds number separately for different angles of attack considered in this work. The results shown in this graph are for a turbulence intensity of 0.25%. As can be seen from this figure, the drag coefficients differ very little for small angles of attack and Reynolds numbers larger than 100,000. Moreover, $C_D\left(Re\right)$ dependencies are not linear functions. The functions that best approximate these calculated coefficient values are third-degree polynomials.

3.2. Turbulence Intensity Effect on Lift and Drag Coefficients

As was reported in Section 3.1, relatively small differences in the Reynolds number result in significant variations in the aerodynamic loads of the airfoil. The effect of turbulence intensity on the aerodynamic characteristics (aerodynamic force coefficients as a function of the angle of attack) of the NACA 0018 airfoil, is much smaller compared to the Reynolds number effect. Nevertheless, the influence of turbulence intensity is particularly important with the lowest value of Re, investigated in this paper, equal to 50,000. Therefore, although the effect of the Reynolds number has been discussed in detail above, this subsection also takes this effect into account. The primary purpose of this subsection is to analyze the effects of turbulence intensity on airfoil performance.

3.2.1. Lift Force Coefficient

Figure 12 shows the characteristics of the lift coefficient, $C_L$, as a function of angle of attack, $\alpha$. These results are compared for the different turbulence intensities considered in this work. Additionally, these results are presented separately for each Reynolds number used in these investigations. The same scale is used in all the graphs in this figure so that the obtained $C_L$ results can be better compared.

It can be seen comparing Figure 6 and Figure 12 that the influence of the turbulence intensity on the lift coefficient is much smaller than the Reynolds number effect. All the $C_L$ characteristics show that the increase in turbulence intensity causes the lift coefficient to decrease in the first region and increase in the second region. In the case of Reynolds numbers larger than 50,000, the significant differences in the $C_L$ values are observed approximately in the middle of the analyzed range of angles of attack. Moreover, as the Reynolds number increases, the differences in $C_L$ results decrease significantly in the second region.

As mentioned in Section 3.1, all simulations carried out in these investigations, regardless of the Reynolds number and the turbulent intensity, showed a zero-lift coefficient at a zero angle of attack.

As discussed above, the turbulence intensity affects the values of the lift coefficients differently in both the first and second regions. The graphs shown in Figure 12 show a common trend for all the calculated series. The turbulence intensity significantly influences the $d{C}_L/d\alpha$ derivative in both regions. This subsection also presents a qualitative analysis of the $d{C}_L/d\alpha$ derivatives for both of these regions. The derivatives of $d{C}_L/d\alpha$ as a function of Reynolds number and for different turbulence intensity values for the first and second regions are shown in Figure 13a,b, respectively. Additionally, Figure 13c also illustrates the transition angle, ${\alpha }_t$, as a function of Reynolds number. The method of determining this angle for this work is also discussed in the previous subsection and Figure 9a. The graphs shown in Figure 13 quantitatively prove three things:

• 1.. An increase in TI causes a decrease in the $d{C}_L/d\alpha$ derivative in the first area, an increase in the $d{C}_L/d\alpha$ derivative in the second area, and an increase in ${\alpha }_t$.

• 2.. An increase in Re causes a decrease in the $d{C}_L/d\alpha$ derivative in the first region and an almost linear increase in the ${\alpha }_t$ angle.

• 3.. Figure 13c does not show a clear trend of the $d{C}_L/d\alpha$ derivative as a function of Re in the second area, both downward and upward. A similar trend of this derivative for the second region is also seen in the experimental studies taken from the literature [79], and is shown in Figure 8.

3.2.2. Drag Force

The results of the drag coefficient, $C_D$, seem to be even less sensitive to the turbulence intensity than the lift coefficients, except for the case of Re = 50,000 (Figure 14). At this Reynolds number, the drag coefficient in the entire range of the angles of attack significantly depends on both the intensity of turbulence and the angle of attack. However, in the case of Re = 200,000, the differences in the drag coefficient results are minimal. In general, an increase in the angle of attack causes an increase in drag. However, as shown in Figure 14, the drag increase is very large at Re = 50,000. In order to compare the results, all the graphs shown in this figure have the same scale.

As mentioned above, the drag coefficient is dependent not only on the Reynolds number and the turbulence intensity but also on the angle of attack. For Reynolds numbers in the range of 100,000 to 200,000, this effect is relatively small. Therefore, Figure 15b only shows the values of the drag coefficient as a function of turbulence intensity for an exemplary Reynolds number of 150,000 because, in the case of Re numbers of 100,000 and 200,000, both the values of the coefficients and their decreases as a function of TI are very similar. Only the results of the $C_D$ coefficient calculated for Re = 50,000 differ significantly from the others. Therefore, these results are shown in a separate figure (Figure 15a). The results shown in Figure 15 show a non-linear relationship of the drag coefficient with the turbulence intensity. Although it is worth noting that it is non-linear, mostly close to TI = 0.1%. In the rest of the TI values taken into account in this paper, an almost constant decrease in the $C_D$ ratio can be observed.

3.3. Static Pressure

The characteristics of the aerodynamic forces discussed in the subsections above result from the pressure distributions around the airfoil. The phenomena occurring in the boundary layer have a significant influence on these distributions. This subsection discusses the effect of the Reynolds number and the turbulence intensity on the static pressure characteristics. In order to better visualize the cause of the occurrence of two regions visible in the lift force coefficient characteristics, the authors of this paper decided to analyze the pressure distributions for both sides of the airfoil separately. Therefore, Figure 16 and Figure 17 show the static pressure coefficients for the suction and pressure sides, respectively. These charts compare the pressure coefficient distributions for different angles of attack and fixed Reynolds numbers. The distributions of the $C_P$ coefficient for the angles of attack closest to the transition angles of attack are shown in black on each of these graphs. The shades of blue represent the angles of attack lower than the transition angle of attack, while the other colors represent the higher angles. It can be seen from Figure 17 that for angles of attack larger than ${\alpha }_t$, there is no longer a laminar-separation bubble on the pressure surface of the airfoil. Moreover, above this angle of attack, the distributions $C_P$ on this side of the profile are already much less dependent on the angle of attack. Figure 16 and Figure 17 show that the largest differences in the $C_P$ coefficients, both on the suction and pressure side of the profile, occur for the smallest Reynolds number of 50,000 analyzed in this study. This explains such significant differences in both the drag coefficient (Figure 10) and the lift coefficient (Figure 6).

Figure 16 and Figure 17 show the effect of the angle of attack for given Reynolds number values. As a result, a qualitative difference was found in the characteristics of the $C_P\left(x/c\right)$ coefficient for the first and second regions (the method of determining the regions is discussed in Section 3.1 of this paper). Figure 18 compares the static pressure coefficients for two angles of attack, one ($\alpha =2°$) for the first region and one ($\alpha =10°$) for the second region. The top two plots shown in this figure compare the $C_P$ relationships for one turbulence intensity value of 0.25%, and the bottom two for one Reynolds number value of 150,000. These characteristics clearly show that the Reynolds number has a much greater effect on the pressure distributions than the turbulence intensity. The effect of both quantities, Reynolds number and turbulence intensity, is much greater on the suction side of the airfoil.

4. Conclusions

This paper discusses the effect of the Reynolds number and the turbulence intensity on the performance of the NACA 0018 profile under the low Reynolds numbers. The research was carried out for the range of Reynolds numbers from 50,000 to 200,000 and for the range of turbulence intensity on the airfoil from 0.01% to 0.5%. Moreover, the investigations were realized for the range of angles of attack from 0 to 10 degrees. Numerical analyzes were performed utilizing the URANS approach with the uncalibrated transition turbulence model $\gamma -R{e}_{\theta }$. Based on the numerical investigations of the airfoil, the following conclusions are drawn:

• One of the key achievements of this paper is the analysis of the aerodynamic characteristics of the NACA 0018 airfoil in relation to the turbulence intensity. The decay rate of turbulence intensity strongly depends on the turbulence intensity at the inlet. The higher the turbulence intensity at the inlet, the higher the decay rate. For the lowest turbulence intensity analyzed in this work, equal to 0.1%, the decay rate is almost constant. In the vicinity of the nose of the profile, the turbulence intensity is very low regardless of the inlet’s turbulence intensity. The results shown in Figure 4b,c have practical applications when using the mesh described in the paper [6]. The data shown in these graphs can be used to interpolate the turbulence intensity at the inlet in order to obtain a specific value of the turbulence intensity on the airfoil.

• All calculations carried out in this paper, regardless of the Reynolds number and the turbulent intensity, showed a zero-lift coefficient at a zero angle of attack.

• The presence of a wide laminar-separation bubble in the boundary layer results in the occurrence of two aerodynamic derivatives of the lift coefficient in the range of angles of attack before stall. The increase in Reynolds number primarily causes a linear increase in the first region. On the other hand, the aerodynamic derivative in the second region is almost independent of the Reynolds number.

• At a Reynolds number of 50,000, a faster increase in drag is observed starting from an angle of attack of 6 deg. For higher Reynolds numbers, the increase in drag coefficient as a function of the angle of attack is much smoother. The drag coefficients differ very little for low angles of attack and Reynolds numbers higher than 100,000.

• Generally, the effect of the turbulence intensity on the lift coefficient is much lower than the Reynolds number effect. The increase in turbulence intensity causes the lift coefficient to decrease in the first region and increase in the second region. Differences in lift coefficients for high angles of attack decrease with increasing the Reynolds number and the turbulence intensity.

• The increase in turbulence intensity causes a decrease in the aerodynamic derivative $d{C}_L/d\alpha$ in the first region and its increase in the second region. Increasing turbulence intensity also increases the transition angle of attack. This trend is visible for all the Reynolds numbers studied in this work.

• For angles of attack higher than the transition angle of attack, there is no longer a laminar-separation bubble on the pressure surface of the airfoil.

• Drag coefficients are less sensitive to the turbulence intensity than the lift coefficients, except for the lowest investigated Reynolds number of 50,000. In general, an increase in the angle of attack causes an increase in the drag coefficient.

• The Reynolds number has a much greater effect on the pressure distributions than the turbulence intensity. The effect of Reynolds number and turbulence intensity is much greater on the suction side of the profile.

• The performed numerical analyzes showed clear differences in the characteristics of the drag coefficient for a Reynolds number of 50,000. For the other considered Reynolds numbers in the range of 100,000 to 200,000, the changes in the drag coefficients are rather linear as a function of both the turbulence intensity and the Reynolds number.

Author Contributions

Conceptualization, K.R.; methodology, J.M. and K.R.; software, J.M. and K.R.; validation J.M. and K.R.; numerical analysis, J.M.; writing—Original draft preparation, K.R.; writing—review and editing, J.M. and K.R.; visualization, J.M. and K.R., supervision, K.R. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

$c$	chord length
$C_D$	drag coefficient
$C_L$	lift coefficient
$C_P$	static pressure coefficient
$d{C}_L/d\alpha$	rate of change of lift coefficient with angle of attack
$l$	turbulence length scale
$L$	length from the airfoil’s trailing edge to the domain’s outlet boundary
$Re$	chord-based Reynolds number
$TI$	turbulence intensity on the airfoil
$T{I}_0$	inlet turbulence intensity
$V_0$	free stream wind velocity [m/s]
$\alpha$	angle of attack
${\alpha }_t$	transition angle of attack (the angle between the first and the second region)
${\delta }_{99}$	boundary-layer thickness

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Table

Figure 1. Silhouette of the H-type Darrieus wind turbine.

Figure 2. Computational domain (a) and mesh around the airfoil (b).

Figure 3. The computing domain and the distance from the trailing edge to the outlet, L (a). The influence of the domain length on the drag coefficient (b) and the lift coefficient (c).

Figure 4. Segment definition for calculating the dependence [Forumla omitted. See PDF.], where [Forumla omitted. See PDF.]—segment start point and TI—segment endpoint (a), turbulence intensity as a function of distance from the airfoil nose (b,c), coefficients of aerodynamic forces as a function of time for different values of turbulence intensity (d–g).

Figure 5. Aerodynamic characteristics of the NACA 0018 profile: lift force coefficient as a function of the angle of attack (a) and drag coefficient (b). Comparison of CFD results with the experiment [24].

Figure 6. Lift coefficient vs. angle of attack. Reynolds number effect.

Figure 7. Two regions of linear growth of lift coefficient curve.

Figure 8. [Forumla omitted. See PDF.] derivatives for two regions. The comparison between present results and data available in the literature [79].

Figure 9. Lift coefficient vs. angle of attack for [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.]. An approach for determining the transition angle, [Forumla omitted. See PDF.] (a). The transition angle vs. Reynolds number for [Forumla omitted. See PDF.] (b).

Figure 10. Drag coefficient vs. angle of attack. Reynolds number effect.

Figure 11. Drag coefficient as a function of Reynolds number.

Figure 12. Lift coefficient vs. angle of attack. Turbulence intensity effect.

Figure 13. The aerodynamic derivatives [Forumla omitted. See PDF.] as a function of Reynolds number for the first region (a) and for the second region (b). Transition angle of attack as a function of Reynolds number (c).

Figure 14. Drag coefficient vs. angle of attack. Turbulence intensity effect.

Figure 15. Effect of turbulence intensity on drag coefficient: (a) for Reynolds number of 50,000 and (b) for Reynolds number of 150,000.

Figure 16. Static pressure coefficient on the suction side at different Reynolds numbers; angle of attack effect.

Figure 17. Static pressure coefficient on pressure side at different Reynolds numbers; angle of attack effect.

Figure 18. Static pressure coefficient on the pressure side (solid lines) and the suction side (dashed lines). Reynolds number effect (upper plots); turbulence intensity effect (lower plots).

Table 1

Turbulence intensity results for mesh sensitivity studies.

	Number of Mesh Elements	TI on the Airfoil [%] 1
Mesh		AoA = 4°	AoA = 10°
Extra coarse	600,000	−0.00041	−0.00041
Coarse	700,000	−0.00016	−0.00016
Medium	800,000	0.24916	0.25065
Fine	900,000	0.00011	0.00011
Extra fine	1,000,000	0.00020	0.00020

¹ Only absolute values for the “Medium” case are presented. For others, differences from the “Medium” case.

References

1. Lin, J.M.; Pauley, L.L. Low-Reynolds-number separation on an airfoil. AIAA J.; 1996; 34, pp. 1570-1577. [DOI: https://dx.doi.org/10.2514/3.13273]

2. Bangga, G.; Dessoky, A.; Wu, Z.; Rogowski, K.; Hansen, M.O.L. Accuracy and consistency of CFD and engineering models for simulating vertical axis wind turbine loads. Energy; 2020; 206, 118087. [DOI: https://dx.doi.org/10.1016/j.energy.2020.118087]

3. Ramírez, J.M.; Saravia, M. Assessment of Reynolds-averaged Navier–Stokes method for modeling the startup regime of a Darrieus rotor. Phys. Fluids; 2021; 33, 037125. [DOI: https://dx.doi.org/10.1063/5.0045083]

4. Branlard, E.; Brownstein, I.; Strom, B.; Jonkman, J.; Dana, S.; Baring-Gould, E.I. A multipurpose lifting-line flow solver for arbitrary wind energy concepts. Wind Energ. Sci.; 2022; 7, pp. 455-467. [DOI: https://dx.doi.org/10.5194/wes-7-455-2022]

5. Silva, J.E.; Danao, L.A.M. Varying VAWT Cluster Configuration and the Effect on Individual Rotor and Overall Cluster Performance. Energies; 2021; 14, 1567. [DOI: https://dx.doi.org/10.3390/en14061567]

6. Rogowski, K.; Królak, G.; Bangga, G. Numerical Study on the Aerodynamic Characteristics of the NACA 0018 Airfoil at Low Reynolds Number for Darrieus Wind Turbines Using the Transition SST Model. Processes; 2021; 9, 477. [DOI: https://dx.doi.org/10.3390/pr9030477]

7. Paraschivoiu, I. Wind Turbine Design: With Emphasis on Darrieus Concept; Polytechnic International Press: Montreal, QC, Canada, 2002.

8. Martinez, J.; Bernabini, L.; Probst, O.; Rodriguez, C. An improved BEM model for the power curve prediction of stall-regulated wind turbines. Wind Energy; 2005; 8, pp. 385-402. [DOI: https://dx.doi.org/10.1002/we.147]

9. Hansen, M.O.L. Aerodynamics of Wind Turbines; Earthscan Publications Ltd.: London, UK, 2008.

10. Zhang, Z.; Nielsen, S.R.K.; Blaabjerg, F.; Zhou, D. Dynamics and Control of Lateral Tower Vibrations in Offshore Wind Turbines by Means of Active Generator Torque. Energies; 2014; 7, pp. 7746-7772. [DOI: https://dx.doi.org/10.3390/en7117746]

11. Fernandez-Gamiz, U.; Zulueta, E.; Boyano, A.; Ansoategui, I.; Uriarte, I. Five Megawatt Wind Turbine Power Output Improvements by Passive Flow Control Devices. Energies; 2017; 10, 742. [DOI: https://dx.doi.org/10.3390/en10060742]

12. Żugaj, M.; Bibik, P.; Jacewicz, M. UAV aircraft model for control system failures analysis. J. Theor. Appl. Mech.; 2016; 54, pp. 1405-1415. [DOI: https://dx.doi.org/10.15632/jtam-pl.54.4.1405]

13. Jimenez, P.; Lichota, P.; Agudelo, D.; Rogowski, K. Experimental Validation of Total Energy Control System for UAVs. Energies; 2020; 13, 14. [DOI: https://dx.doi.org/10.3390/en13010014]

14. Sheng, H.; Zhang, C.; Xiang, Y. Mathematical Modeling and Stability Analysis of Tiltrotor Aircraft. Drones; 2022; 6, 92. [DOI: https://dx.doi.org/10.3390/drones6040092]

15. El-Fahham, I.; Abdelshahid, G.; Mokhiamar, O. Pitch Angle Modulation of the Horizontal and Vertical Axes Wind Turbine Using Fuzzy Logic Control. Processes; 2021; 9, 1337. [DOI: https://dx.doi.org/10.3390/pr9081337]

16. Sakamoto, L.; Fukui, T.; Morinishi, K. Blade Dimension Optimization and Performance Analysis of the 2-D Ugrinsky Wind Turbine. Energies; 2022; 15, 2478. [DOI: https://dx.doi.org/10.3390/en15072478]

17. Ghiasi, P.; Najafi, G.; Ghobadian, B.; Jafari, A.; Mazlan, M. Analytical Study of the Impact of Solidity, Chord Length, Number of Blades, Aspect Ratio and Airfoil Type on H-Rotor Darrieus Wind Turbine Performance at Low Reynolds Number. Sustainability; 2022; 14, 2623. [DOI: https://dx.doi.org/10.3390/su14052623]

18. Ferreira, C.S.; Madsen, H.A.; Barone, M.; Roscher, B.; Deglaire, P.; Arduin, I. Comparison of aerodynamic models for Vertical Axis Wind Turbines. J. Phys. Conf. Ser.; 2014; 524, 012125. [DOI: https://dx.doi.org/10.1088/1742-6596/524/1/012125]

19. Castelein, D.; Ragni, D.; Tescione, G.; Ferreira, C.J.S.; Gaunaa, M. Creating a benchmark of Vertical Axis Wind Turbines in Dynamic Stall for validating numerical models. Proceedings of the 33rd Wind Energy Symposium, AIAA SciTech, AIAA 2015-0723; Kissimmee, FL, USA, 5–9 January 2015; [DOI: https://dx.doi.org/10.2514/6.2015-0723]

20. Rogowski, K.; Hansen, M.O.L.; Bangga, G. Performance Analysis of a H-Darrieus Wind Turbine for a Series of 4-Digit NACA Airfoils. Energies; 2020; 13, 3196. [DOI: https://dx.doi.org/10.3390/en13123196]

21. Lichota, P. Multi-Axis Inputs for Identification of a Reconfigurable Fixed-Wing UAV. Aerospace; 2020; 7, 113. [DOI: https://dx.doi.org/10.3390/aerospace7080113]

22. Lichota, P.; Norena, D.A. A Priori Model Inclusion in the Multisine Maneuver Design. Proceedings of the 17th International Carpathian Control Conference; High Tatras, Slovakia, 29 May–1 June 2016; Petras, I.; Podlubny, I.; Kacur, J. IEEE: Tatranská Lomnica, Slovakia, 2016; pp. 440-445. [DOI: https://dx.doi.org/10.1109/CarpathianCC.2016.7501138]

23. ANSYS, Inc. ANSYS Fluent Theory Guide, Release 19.1.; ANSYS, Inc.: Canonsburg, PA, USA, 2019.

24. Timmer, W.A. Two-Dimensional Low-Reynolds Number Wind Tunnel Results for Airfoil NACA 0018. Wind. Eng.; 2008; 32, pp. 525-537. [DOI: https://dx.doi.org/10.1260/030952408787548848]

25. Dul, F.; Lichota, P.; Rusowicz, A. Generalized Linear Quadratic Control for a Full Tracking Problem in Aviation. Sensors; 2020; 20, 2955. [DOI: https://dx.doi.org/10.3390/s20102955]

26. Ren, G.; Liu, J.; Wan, J.; Li, F.; Guo, Y.; Yu, D. The analysis of turbulence intensity based on wind speed data in onshore wind farms. Renew. Energy; 2018; 123, pp. 756-766. [DOI: https://dx.doi.org/10.1016/j.renene.2018.02.080]

27. Laneville, A.; Vittecoq, P. Dynamic stall: The case of the vertical axis wind turbine. J. Sol. Energy Eng.; 1986; 108, pp. 140-145. [DOI: https://dx.doi.org/10.1115/1.3268081]

28. Rogowski, K.; Maroński, R.; Hansen, M.O.L. Steady and unsteady analysis of NACA 0018 airfoil in vertical-axis wind turbine. J. Theor. Appl. Mech.; 2018; 51, pp. 203-212. [DOI: https://dx.doi.org/10.15632/jtam-pl.56.1.203]

29. Yang, Y.; Guo, Z.; Zhang, Y.; Jinyama, H.; Li, Q. Numerical Investigation of the Tip Vortex of a Straight-Bladed Vertical Axis Wind Turbine with Double-Blades. Energies; 2017; 10, 1721. [DOI: https://dx.doi.org/10.3390/en10111721]

30. Mohan Kumar, P.; Sivalingam, K.; Lim, T.-C.; Ramakrishna, S.; Wei, H. Strategies for Enhancing the Low Wind Speed Performance of H-Darrieus Wind Turbine—Part 1. Clean Technol.; 2019; 1, pp. 185-204. [DOI: https://dx.doi.org/10.3390/cleantechnol1010013]

31. Rogowski, K. CFD Computation of the H-Darrieus Wind Turbine—The Impact of the Rotating Shaft on the Rotor Performance. Energies; 2019; 12, 2506. [DOI: https://dx.doi.org/10.3390/en12132506]

32. Rogowski, K.; Hansen, M.O.L.; Maroński, R.; Lichota, P. Scale Adaptive Simulation Model for the Darrieus Wind Turbine. J. Phys. Conf. Ser.; 2016; 753, 022050. [DOI: https://dx.doi.org/10.1088/1742-6596/753/2/022050]

33. Alaimo, A.; Esposito, A.; Messineo, A.; Orlando, C.; Tumino, D. 3D CFD Analysis of a Vertical Axis Wind Turbine. Energies; 2015; 8, pp. 3013-3033. [DOI: https://dx.doi.org/10.3390/en8043013]

34. Jacobs, N.E.; Ward, K.E.; Pinkerton, R.M. The Characteristics of 78 Related Airfoil Sections from Tests in the Variable-Density Wind Tunnel. NACA Report 460; Langley Memorial Aeronautical Laboratory: Hampton, VA, USA, 1933; ISBN NACA-460 Available online: https://ntrs.nasa.gov/api/citations/19930091108/downloads/19930091108.pdf (accessed on 15 April 2022).

35. Anderson, J.D., Jr. Fundamentals of Aerodynamics; 5th ed. McGraw-Hill: New York, NY, USA, 2011.

36. Moran, J. An Introduction to Theoretical and Computational Aerodynamics; 1st ed. John Wiley & Sons: New York, NY, USA, 1984.

37. Królak, G. Numerical Analysis of Aerodynamic Characteristics of the NACA 0018 Airfoil Used in the Darrieus wind Turbines. Master’s Thesis; The Faculty of Power and Aeronautical Engineering, Warsaw University of Technology: Warsaw, Poland, 2020.

38. Ferreira, C.S.; Barone, M.; Zanon, A.; Giannattasio, P. Airfoil optimization for stall regulated vertical axis wind turbines. Proceedings of the AIAA SciTech—33rd Wind Energy Symposium; Kissimmee, FL, USA, 5–9 January 2015; pp. 1-16. [DOI: https://dx.doi.org/10.2514/6.2015-0722]

39. Ashwill, T.D. Measured Data for the Sandia 34-Meter Vertical Axis Wind Turbine; SAND91-2228 Sandia National Laboratories: Albuquerque, NM, USA, 1992.

40. Shires, A. Development and Evaluation of an Aerodynamic Model for a Novel Vertical Axis Wind Turbine Concept. Energies; 2013; 6, pp. 2501-2520. [DOI: https://dx.doi.org/10.3390/en6052501]

41. Jacobs, E.N.; Sherman, A. Airfoil Section Characteristics as Affected by Variations of the Reynolds Number; Technical Report National Advisory Committee for Aeronautics: Langley Field, VA, USA, 1937; Available online: https://ntrs.nasa.gov/api/citations/19930091662/downloads/19930091662.pdf (accessed on 15 April 2022).

42. Abbott, H.; Doenhoff, E.; Lous, S.; Stivers, J. National Advisory Commitee for Aeronautics Report No. 824: Summary of Airfoil Data. 1945; Available online: https://ntrs.nasa.gov/api/citations/19930090976/downloads/19930090976.pdf (accessed on 15 April 2022).

43. Sheldahl, R.E.; Klimas, P.C. Aerodynamic Characteristics of Seven Symmetrical Airfoil Sections through 180-Degree Angle of Attack for Use in Aerodynamic Analysis of Vertical Axis Wind Turbines; Technical Report Sandia National Laboratories: Albuquerque, NM, USA, 1981; [DOI: https://dx.doi.org/10.2172/6548367]

44. Bastankhah, M.; Porté-Agel, F. A New Miniature Wind Turbine for Wind Tunnel Experiments. Part I: Design and Performance. Energies; 2017; 10, 908. [DOI: https://dx.doi.org/10.3390/en10070908]

45. Hezaveh, S.H.; Bou-Zeid, E.; Lohry, M.W.; Martinelli, L. Simulation and wake analysis of a single vertical axis wind turbine. Wind Energy; 2017; 20, pp. 713-730. [DOI: https://dx.doi.org/10.1002/we.2056]

46. Bedon, G.; Schmidt Paulsen, U.; Aagaard Madsen, H.; Belloni, F.; Raciti Castelli, M.; Benini, E. Computational assessment of the DeepWind aerodynamic performance with different blade and airfoil configurations. Appl. Energy; 2017; 185, pp. 1100-1108. [DOI: https://dx.doi.org/10.1016/j.apenergy.2015.10.038]

47. Rossander, M.; Dyachuk, E.; Apelfröjd, S.; Trolin, K.; Goude, A.; Bernhoff, H.; Eriksson, S. Evaluation of a Blade Force Measurement System for a Vertical Axis Wind Turbine Using Load Cells. Energies; 2015; 8, pp. 5973-5996. [DOI: https://dx.doi.org/10.3390/en8065973]

48. Shires, A.; Kourkoulis, V. Application of Circulation Controlled Blades for Vertical Axis Wind Turbines. Energies; 2013; 6, pp. 3744-3763. [DOI: https://dx.doi.org/10.3390/en6083744]

49. Nakano, T.; Fujisawa, N.; Oguma, Y.; Takagi, Y.; Lee, S. Experimental study on flow and noise characteristics of NACA0018 airfoil. J. Wind. Eng. Ind. Aerodyn.; 2007; 95, pp. 511-531. [DOI: https://dx.doi.org/10.1016/j.jweia.2006.11.002]

50. Bianchini, A.; Balduzzi, F.; Rainbird, J.M.; Peiro, J.; Graham, J.M.R.; Ferrara, G.; Ferrari, L. An Experimental and Numerical Assessment of Airfoil Polars for Use in Darrieus Wind Turbines—Part I: Flow Curvature Effects. J. Eng. Gas Turbine Power; 2016; 138, 032602. [DOI: https://dx.doi.org/10.1115/1.4031269]

51. Bianchini, A.; Balduzzi, F.; Rainbird, J.M.; Peiro, J.; Graham, J.M.R.; Ferrara, G.; Ferrari, L. An Experimental and Numerical Assessment of Airfoil Polars for Use in Darrieus Wind Turbines—Part II: Post-stall Data Extrapolation Methods. J. Eng. Gas. Turbine Power; 2016; 138, 032603. [DOI: https://dx.doi.org/10.1115/1.4031270]

52. Istvan, M. Turbulence intensity effects on laminar separation bubbles formed over an airfoil. AIAA J.; 2018; 56, pp. 335-1347. [DOI: https://dx.doi.org/10.2514/1.J056453]

53. Istvan, M.S.; Yarusevych, S. Effects of free-stream turbulence intensity on transition in a laminar separation bubble formed over an airfoil. Exp. Fluids; 2018; 59, 52. [DOI: https://dx.doi.org/10.1007/s00348-018-2511-6]

54. Claessens, M.C. The Design and Testing of Airfoils in Small Vertical Axis Wind Turbines. Master’s Thesis; Delft University of Technology: Delft, The Netherlands, 2006; Available online: https://repository.tudelft.nl/islandora/object/uuid%3Afe4a7ae3-103e-40f3-a1f1-ae6d910d8c71 (accessed on 10 April 2022).

55. Du, L. Numerical and Experimental Investigations of Darrieus Wind Turbine Start-up and Operation. Ph.D. Thesis; Durham University: Durham, UK, 2016; Available online: http://etheses.dur.ac.uk/11384/ (accessed on 10 April 2022).

56. Müller-Vahl, H.F.; Strangfeld, C.; Nayeri, C.; Paschereit, C.O.; Greenblatt, D. Control of Thick Airfoil, Deep Dynamic Stall Using Steady Blowing. AIAA J.; 2015; 53, pp. 277-295. [DOI: https://dx.doi.org/10.2514/1.J053090]

57. Melani, P.F.; Balduzzi, F.; Ferrara, G.; Bianchini, A. An Annotated Database of Low Reynolds Aerodynamic Coefficients for the NACA0018 Airfoil. AIP Conf. Proc.; 2019; 2191, 020110. [DOI: https://dx.doi.org/10.1063/1.5138843]

58. Belabes, B.; Paraschivoiu, M. Numerical study of the effect of turbulence intensity on VAWT performance. Energy; 2021; 233, 121139. [DOI: https://dx.doi.org/10.1016/j.energy.2021.121139]

59. Tangermann, E.; Klein, M. Numerical simulation of laminar separation on a NACA0018 airfoil in freestream turbulence. Proceedings of the AIAA Scitech 2020 Forum; Orlando, FL, USA, 6–10 January 2020; pp. 1-10. [DOI: https://dx.doi.org/10.2514/6.2020-2064]

60. Marxen, O.; Rist, U.; Wagner, S. Effect of spanwise-modulated disturbances on transition in a separated boundary layer. AIAA J.; 2004; 42, pp. 937-944. [DOI: https://dx.doi.org/10.2514/1.565]

61. Galbraith, M.C.; Visbal, M.R. Implicit Large Eddy Simulation of Low-Reynolds-Number Transitional Flow Past the SD7003 Airfoil. Proceedings of the 40th Fluid Dynamics Conference and Exhibit; Chicago, IL, USA, 28 June–1 July 2010; [DOI: https://dx.doi.org/10.2514/6.2010-4737]

62. Zhiyin, Y. Large-eddy simulation: Past, present and the future. Chin. J. Aeronaut.; 2015; 28, pp. 11-24. [DOI: https://dx.doi.org/10.1016/j.cja.2014.12.007]

63. Kim, J.-H.; Ahn, J. Large Eddy Simulation of Leakage Flow in a Stepped Labyrinth Seal. Processes; 2021; 9, 2179. [DOI: https://dx.doi.org/10.3390/pr9122179]

64. Huang, Z.; Shi, G.; Liu, X.; Wen, H. Effect of Flow Rate on Turbulence Dissipation Rate Distribution in a Multiphase Pump. Processes; 2021; 9, 886. [DOI: https://dx.doi.org/10.3390/pr9050886]

65. Walters, D.K.; Cokljat, D. A three-equation eddy-viscosity model for Reynolds-averaged Navier-Stokes simulations of transitional flow. J. Fluids Eng.; 2008; 130, 121401. [DOI: https://dx.doi.org/10.1115/1.2979230]

66. Dick, E.; Kubacki, S. Transition Models for Turbomachinery Boundary Layer Flows: A Review. Int. J. Turbomach. Propuls. Power; 2017; 2, 4. [DOI: https://dx.doi.org/10.3390/ijtpp2020004]

67. Menter, F.R.; Langtry, R.B.; Likki, S.R.; Suzen, Y.B.; Huang, P.G.; Völker, S. A Correlation-Based Transition Model Using Local Variables: Part I—Model Formulation. J. Turbomach.; 2006; 128, pp. 413-422. [DOI: https://dx.doi.org/10.1115/1.2184352]

68. Turner, C. Laminar Kinetic Energy Modelling for Improved Laminar-Turbulent Transition Prediction. Ph.D. Thesis; School of Mechanical, Aerospace and Civil Engineering, The University of Manchester: Manchester, UK, 2012; Available online: https://www.proquest.com/pagepdf/1774233482 (accessed on 10 April 2022).

69. Suluksna, K.; Juntasaro, E. Assessment of intermittency transport equations for modeling transition in boundary layers subjected to freestream turbulence. Int. J. Heat Fluid Flow; 2008; 29, pp. 48-61. [DOI: https://dx.doi.org/10.1016/j.ijheatfluidflow.2007.08.003]

70. Choudhry, A.; Arjomandi, M.; Kelso, R. A study of long separation bubble on thick airfoils and its consequent effects. Int. J. Heat Fluid Flow; 2015; 52, pp. 84-96. [DOI: https://dx.doi.org/10.1016/j.ijheatfluidflow.2014.12.001]

71. Aftab, S.M.A.; Mohd Rafie, A.S.; Razak, N.A.; Ahmad, K.A. Turbulence Model Selection for Low Reynolds Number Flows. PLoS ONE; 2016; 11, e0153755. [DOI: https://dx.doi.org/10.1371/journal.pone.0153755] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27104354]

72. Michna, J.; Rogowski, K.; Bangga, G.; Hansen, M.O.L. Accuracy of the Gamma Re-Theta Transition Model for Simulating the DU-91-W2-250 Airfoil at High Reynolds Numbers. Energies; 2021; 14, 8224. [DOI: https://dx.doi.org/10.3390/en14248224]

73. Tescione, G.; Ragni, D.; He, C.; Ferreira, C.J.S.; van Bussel, G.J.W. Near wake flow analysis of a vertical axis wind turbine by stereoscopic particle image velocimetry. Renew. Energy; 2014; 70, pp. 47-61. [DOI: https://dx.doi.org/10.1016/j.renene.2014.02.042]

74. Guo, Z.; Fletcher, D.F.; Haynes, B.S. Implementation of a height function method to alleviate spurious currents in CFD modelling of annular flow in microchannels. Appl. Math. Model.; 2015; 39, pp. 4665-4686. [DOI: https://dx.doi.org/10.1016/j.apm.2015.04.022]

75. Kapsalis, P.C.S.; Voutsinas, S.; Vlachos, N.S. Comparing the effect of three transition models on the CFD predictions of a NACA0012 airfoil aerodynamics. J. Wind. Eng. Ind.; 2016; 157, pp. 158-170. [DOI: https://dx.doi.org/10.1016/j.jweia.2016.07.007]

76. Morgado, J.; Vizinho, R.; Silvestre, M.A.R.; Páscoa, J.C. XFOIL vs. CFD performance predictions for high lift low Reynolds number airfoils. Aerosp. Sci. Technol.; 2016; 52, pp. 207-214. [DOI: https://dx.doi.org/10.1016/j.ast.2016.02.031]

77. Günel, O.; Koç, E.; Yavuz, T. CFD vs. XFOIL of airfoil analysis at low reynolds numbers. In Proceedings of the 2016 IEEE International Conference on Renewable Energy Research and Applications (ICRERA); Birmingham, UK, 20–23 November 2016; pp. 628-632. [DOI: https://dx.doi.org/10.1109/ICRERA.2016.7884411]

78. Cakmakcioglu, S.C.; Sert, I.O.; Tugluk, O.; Sezer-Uzol, N. 2-D and 3-D CFD Investigation of NREL S826 Airfoil at Low Reynolds Numbers. J. Phys. Conf. Ser.; 2014; 524, 012028. [DOI: https://dx.doi.org/10.1088/1742-6596/524/1/012028]

79. Gerakopulos, R.; Boutilier, M.S.H.; Yarusevych, S. Aerodynamic Characterization of a NACA 0018 Airfoil at Low Reynolds Numbers. Proceedings of the AIAA 2010-4629, 40th Fluid Dynamics Conference and Exhibit; Chicago, IL, USA, 28 June–1 July 2010; [DOI: https://dx.doi.org/10.2514/6.2010-4629]