An open database of high-fidelity, multi-Reynolds airfoil polars for wind turbine blade design

This paper presents an open-access dataset that has been developed to provide high-fidelity aerodynamic polars for a wide range of wind turbine airfoils, computed using a consistent computational fluid dynamics (CFD) methodology. The database includes lift, drag, and moment coefficients across multiple Reynolds and Mach numbers. Coefficients are computed with both fully turbulent and free transition boundary layers between positive and negative stall angles. A blend of the free transition and fully turbulent coefficients is also included. Beyond the stall region, state-of-the-art post-stall extrapolation models are used, the calibration parameters of which are derived from a number of high-fidelity calculations at various angles of attack in separated flow. The data are relevant for wind turbine design, modeling, and simulation and reflect representative airfoil performance along the span of modern, large-size offshore rotors. The database includes FFA-W3, FFA-W2, FFA-W1, DU, FX-77, and the recently developed OSO families of airfoils. All simulations are performed using validated numerical methods and span a range of conditions typical of utility-scale offshore turbines. The paper discusses the dataset and its key differences from other open-access data, as well as some general trends that can be noted in the aerodynamic coefficients; such a discussion is made possible by the volume of data produced and is specifically tailored to wind turbine applications.

1 Background and motivation

Wind turbines face a unique set of operating conditions during their lifetime (Veers et al., 2023). These machines mostly operate inside the atmospheric boundary layer (ABL), often in a highly turbulent flow, and are subject to spatial and temporal variations in wind speed, wind direction, and temperature. Recent large-sized rotor blades may even partially operate between the ABL and the free atmosphere during events of strong atmospheric stability, facing inflow conditions unprecedented to date. Moreover, these machines are constantly exposed to atmospheric elements while being expected to operate reliably for many years with minimal levels of maintenance. This makes them prone to performance degradation due to soiling and erosion over time.

These challenges are reflected in the evolution of wind turbine blade shapes over the years (van Bussel, 2022). The core geometry of a wind turbine blade consists of radially stacked two-dimensional airfoil sections. Early wind turbine designs from the 1970s generally adopted pre-World War II standard airfoil shapes designed for aircraft wings from the National Advisory Committee for Aeronautics (NACA) (Abbott et al., 1945). These airfoils mostly featured relatively low thickness-to-chord ratios and narrow optimal operating windows. Such characteristics proved ill-suited to the highly unsteady conditions experienced by wind turbines and posed significant limitations on blade scaling. Moreover, they often suffered severe performance degradation due to soiling or small-scale damage, with wind turbine power production dropping nearly 30 % in some cases (Airfoils, Where the Turbine Meets the Wind, 2025). These considerations motivated researchers to develop airfoil shapes tailored specifically to wind turbine applications.

Development efforts were significantly accelerated by the introduction of panel methods, which have been studied since the early 1960s (Erickson, 1990). In particular, the coupling of inviscid panel methods with viscous boundary layer models in the early 1980s made low-speed aerodynamic analyses viable (Drela and Giles, 1987). During this decade the increase in computational power and the development of inverse design approaches, where designers specified the desired pressure distribution and where the panel method would compute an airfoil shape (Drela, 1989), significantly catalyzed the development of wind turbine airfoils. A hallmark of this evolution is that, since its beginning, it advanced in parallel within academia and wind turbine manufacturers. However, due to industrial confidentiality, commercial airfoil geometries were generally not disclosed, slowing down the validation of computational methods and the development of shared databases for wind turbine design.

From the start, airfoil designers developed families of airfoils, comprising thicker geometries which prioritized structural characteristics and more slender shapes to achieve high aerodynamic efficiency. An overview of the evolution of wind turbine airfoils in terms of aerodynamic efficiency and bending stiffness, which is assumed to be proportional to the sectional moment of inertia around the chord, is provided in Fig. 1. The figure is based on publicly available data gathered with different methods and in different testing conditions and which only partially cover the entire family in some cases. It should therefore be used to infer general trends in wind turbine airfoil design rather than as a comparison of the relative strengths of the airfoil families. One of the first airfoil families to be publicly issued is the S series, released by the National Laboratory of the Rockies (NLR; formerly known as NREL – the National Renewable Energy Laboratory) in the 1990s (Tangler and Somers, 1995). S-series airfoils feature multiple families with varying thickness-to-chord ratios and have been used in multiple commercial wind turbines (Thresher et al., 1994). However, as shown in Fig. 1, these airfoils are still relatively thin for the structural requirements of modern 100+ m long blades and are relatively sensitive to changes in roughness. In the same period (late 1970s), the FX-77 family of airfoils was developed at the University of Stuttgart (Wortmann, 1978). This family featured airfoils with relative thicknesses between 15.3 % and 50 %, with the thicker shapes derived from truncating a 34.3 %-thick design (flatback airfoil). To cater to the increased structural requirements of larger wind turbine blades, airfoil designers started developing thicker airfoil shapes, especially for use in inboard regions of the blade, as reflected by the lower efficiency and much higher sectional stiffness of root airfoils shown in Fig. 1. The University of Delft (the Netherlands) has developed several airfoils for wind turbine applications from 1990 onwards (Timmer and Rooij, 2003). They can be easily identified by their name: DU (Delft University) – W (wind) – followed by the year of development and the thickness of the design. Several airfoil families have also been developed by Risø in Denmark, namely the A1 (Bertagnolio et al., 2001; Dahl and Fuglsang, 1998), B1 (Fuglsang et al., 2004), and C2 families (Bak et al., 2008). These airfoils are relatively recent and have been designed with aero-structural considerations in mind (Bak et al., 2014). Risø airfoils, however, are used commercially and are not available open access. The Aeronautical Research Institute of Sweden developed a series of publicly available wind turbine airfoils for wind turbine applications in 1990 (Andres, 1990). They have become very popular among the academic community due to their use in several research reference turbines, such as the DTU 10 MW (Bak et al., 2013), IEA 10 and 3.4 MW (Bortolotti et al., 2019), IEA 15 MW (IEA Wind, 2020), and IEA 22 MW (Zahle et al., 2024). Moreover, they feature high thickness-to-chord ratios and good performance when soiled. Several airfoil families have been designed more recently by various institutions, such as ECN (Boorsma et al., 2015; Grasso, 2014), CENER (Méndez et al., 2014), and others (Cheng et al., 2014; Hansen, 2018), using a variety of design and optimization methods. These airfoil families are, however, not open access, and the shapes of the airfoils are not available publicly. On the other hand, SANDIA has recently developed a family of airfoils specifically tailored to offshore wind turbines (Karcher et al., 2025), i.e., for 100+ m long blades with high Reynolds (Re) and Mach (Ma) numbers. Indeed, the sheer size of modern wind turbine blades has increased the Reynolds numbers dramatically, especially offshore. Very limited experimental validation data exist at these Reynolds numbers due to experimental limitations, which makes the calibration and validation of simulation tools challenging.

Figure 1 Historical evolution of airfoil families for wind turbine design in terms of design balance between aerodynamic efficiency and structural bending stiffness, assumed to be proportional to the flapwise sectional moment of inertia. The performance envelope represents the range between tip (clean) and root (rough) profiles for each family depending on the availability of public data: NACA 63-415/63(4)421 (Bertagnolio et al., 2001), NREL S801/S815 (Ramsay et al., 1996), FFA-W3 211/301 (Andres, 1990), DU 96-W-180/97-W-300 (Timmer and Rooij, 2003), Risø A1-18/24 (Fuglsang et al., 1999), Risø B1-18/24 (Fuglsang et al., 2003), and OSO 21/30 (Karcher et al., 2025). Percentages refer to relative thickness of the selected root and tip airfoils.

Despite research showing how these rotating blades feature strong three-dimensional flow characteristics, especially near the blade root (Guntur and Sørensen, 2015; Hand et al., 2001), modeling them as a series of bi-dimensional sections allows for several engineering simplifications and is the core approach used in most modeling tools. As such, all modern aero-elastic tools are based on polar data, and therefore accurate airfoil aerodynamic coefficients are essential to obtain reliable results. Because airfoil-based models are used in some capacity along the entire design chain of a wind turbine, accuracy impacts all phases of the process, from initial blade design to load verification and certification. In the absence of solid experimental reference data, the aerodynamic coefficients of a given airfoil can vary greatly based on the method used to compute the coefficients and the assumptions made. For instance, the aerodynamic coefficients used in the design of the DTU 10 MW, the IEA 15 MW, and the IEA 22 MW differ, despite all three designs featuring the same FFA-W3 family of airfoils. These uncertainties also make comparisons between different airfoil families difficult. Due to being developed in different times and using different tools, the aerodynamic coefficients used in the design of the aforementioned airfoil families cannot be directly compared in most cases. Moreover, even for those families that have been extensively tested experimentally, such as the DU series of airfoils (Timmer and Rooij, 2003), the FX-77 (Wortmann, 1978), some of the FFA-W3 airfoils (Fuglsang et al., 1998), or the S-series airfoils (Ohio State University Wind Tunnel Tests | Wind Research, 2026), making direct comparisons is challenging due to the coefficients being computed by different researchers in different wind tunnels, with different surface finishes, inflow turbulence, blockage ratios, and Reynolds and Mach numbers. In fact, the complex physics at play cause the performance of an aerodynamic shape to be strongly influenced by minor variations in external shape or surface finish and changes in inflow conditions.

Past efforts to gather coefficients for multiple airfoil families include the ECN (now TNO) Aerodynamic Table Generator, which utilizes wind tunnel measurements and empirical models to populate aerodynamic databases for industrial design codes. While historically significant, the tool is documented in Dutch and remains a proprietary institutional asset, limiting its availability to the broader research community (Bot, 2001). More recently, the viscous panel method RFoil (Koodly Ravishankara et al., 2025) was also used to provide coefficients for the FFA-W3, OSO, RISO, and DU families of airfoils in the work led by Sandia National Laboratories (Maniaci et al., 2025).

In addition, significant uncertainty is also present in the post-stall region. In fact, aerodynamic coefficients in the attached flow regime are typically obtained using viscous panel methods, computational fluid dynamics (CFD), or experimental measurements at moderate angles of attack. In contrast, post-stall extrapolation methods are widely used to empirically extend the coefficients to a −180 to +180° range. Once again, depending on the method and the assumptions used, estimations can vary greatly. For instance, the maximum drag value within the coefficients used in the NREL 5 MW RWT designs is approximately 1.3 before stall-delay correction, while this value is 1.3 to 1.5 for the DTU 10 MW, IEA 15 MW, and IEA 22 MW, depending on the airfoil thickness. As wind turbine blades have complex geometries and a finite aspect ratio, these values are often set empirically based on experience or previous simulations (Zahle et al., 2024). Indeed, wind tunnel experiments on two-dimensional semi-infinite sections, such as the ones presented in Timmer (2020), indicate maximum drag coefficients above 1.8 for a selection of wind turbine airfoils. Similar values are also recommended in some of the original post-extrapolation models, like that of Viterna and Janetzke (1982), which suggest a value of approximately 1.8 depending on the specific characteristics of the airfoil such as camber, relative thickness, and leading-edge radius. These differences cannot be imputed to Reynolds number alone (Battisti et al., 2020) but rather to the finite aspect ratio of wind turbine blades and their complex geometry, as noted in Burton et al. (2011), where mean drag coefficients measured on standstill blades are discussed. Overall, the values that should be used in engineering models remain highly uncertain. Although horizontal-axis wind turbines do not typically operate in the post-stall region, these conditions can be encountered when parked and can give rise to unsteady phenomena, such as vortex-induced vibration (VIV), which can drive ultimate loading. Even in the absence of unsteady loading on the blade, accurate sectional coefficients at high angles of attack can help reduce modeling errors in engineering tools.

Due to the issues described so far, computing reliable airfoil coefficients can be challenging for wind energy practitioners who wish to build their own blade design, as it requires significant expertise. It is, in fact, not uncommon to see coefficients computed using low-order viscous panel methods in recent research efforts (Gupta et al., 2024), introducing significant uncertainty in drag coefficient estimation and in the near-stall region (Ramanujam et al., 2016).

This work aims to address these issues by computing a dataset of high-Reynolds airfoil coefficients that can serve as a common research basis for modern wind turbine design. One of the key aspects of this effort is the fact that the airfoil polars have been computed using the same modeling choices and inflow conditions for the various families included in the dataset, thus ensuring comparability. Moreover, the coefficients included in this dataset take full advantage of CFD and are computed including compressibility, the effects of which are becoming increasingly relevant for modern wind turbines (Vitulano et al., 2025).

Multiple airfoil families are included in the dataset, allowing for fair comparisons between them. Coefficients are computed using consistent state-of-the-art numerical methods, which are described in the following sections of this study.

The study is organized as follows. Section 2 outlines the characteristics of the dataset, including the airfoils and simulated conditions. Section 3 describes the numerical methods used to compute the dataset, including the validation and verification of such methods. Section 4 presents a selection of results.

2 Wind turbine airfoils

The dataset analyzed herein contains airfoil coefficients with varying Reynolds and Mach numbers of several open-access commonly used wind turbine airfoils. The dataset includes the FFA-W1, FFA-W2, FFA-W3, DU, FX77, and OSO airfoil families. One notable exception amongst open-access airfoils is represented by NREL S-series airfoils (Tangler and Somers, 1995), as they are relatively slender for modern blade designs. In fact, only airfoils with thickness-to-chord ratios greater than 21 % are included in the dataset. The NACA63(3)418 airfoil is included despite featuring a relative thickness of 18 % as it has been used in recent actual wind turbine blades (Braud et al., 2025). The full list of geometries is shown in Fig. 2 and is summarized in Table 1.

Table 1 Airfoils included in the dataset. Free transition: simulations with boundary layer transition (clean blade). Fully turbulent: simulations without boundary layer transition (soiled blade). Blend refers to a 70 % free transition and 30 % fully turbulent interpolation, to be used as a reference for blade design.

Figure 2 Airfoils included in the present database grouped by airfoil family.

The Reynolds (Eq. 1) and Mach (Eq. 2) number are defined as

$$\begin{array}{}\text{(1)}& {\displaystyle }Re={\displaystyle \frac{\mathit{\rho }cU}{\mathit{\mu }}},\text{(2)}& {\displaystyle }Ma={\displaystyle \frac{U}{\sqrt{\mathit{\gamma }RT}}},\end{array}$$

where ρ is the air density, μ is the dynamic viscosity of the fluid, γ is the ratio of heat capacity at constant pressure and at constant volume, R is the gas constant, and T is the temperature of the fluid. Modern wind turbine blades, particularly those of offshore machines, are increasingly large. As such, chord-based Reynolds numbers can exceed 15×10⁶ (Karcher et al., 2025). The mean values assumed in this study and their measurement units are summarized in Table 2, which also reports the reference conditions in terms of pressure P and turbulence intensity ti considered for the simulations. A turbulence intensity of 0.35 % is considered to match the conditions where the computational model is tuned and validated, as explained in detail in Sect. 3.1. In addition, this value is in line with the state-of-the-art. For instance, it is not uncommon to see aerodynamic coefficients used in blade design exercise being computed with N_crit values in the range of 7 to 9 (Zahle et al., 2024), which correspond to low levels of turbulence intensity, as shown in Eq. (6) discussed later on in this study.

It must be noted that the inflow turbulence wind turbine airfoils experience during operation varies widely due to many factors. Depending on the turbulence level in the incoming wind field, the actual turbulence on the blade section scales based on the relative velocity, which depends on the local tip speed ratio of each section. Especially in the case of outboard blade sections, this can lower the turbulence intensity level significantly. The many factors at play make choosing a single value challenging. However, providing users of the dataset with coefficients computed with low inflow turbulence and a fully turbulent boundary layer allows for interpolation between the two, in order to fine-tune the coefficients to the application and the level of desired conservativeness, if used for blade design applications.

All the simulations performed in this study model the fluid as compressible. This is a key distinction compared to the aerodynamic coefficients used in common academic reference wind turbines (Bortolotti et al., 2019; IEA Wind, 2020; Zahle et al., 2024), where the flow is modeled as incompressible. This simplification is becoming increasingly limiting, as recent studies have shown how, due to high tip speeds and local flow acceleration, transonic flow is possible in certain conditions near the blade tip (De Tavernier and Terzi, 2022). This condition cannot be appropriately modeled if the flow is considered incompressible. On the other hand, Vitulano et al. (2025) have shown that unsteady Reynolds-averaged Navier–Stokes (URANS) simulations are a viable way of predicting airfoil performance in transonic conditions. In addition, compressibility is known to have an effect on the aerodynamic coefficients, particularly on the slope of the lift coefficient and the minimum drag coefficient in the pre-stall region, which tends to increase as the Mach number increases (Gudmundsson, 2014).

Table 2 Parameters used in the current study.

The second key improvement of this dataset is the fact that multiple Reynolds numbers are considered. In fact, the Reynolds number can vary greatly during operation in actual wind turbines, as the relative inflow velocity changes due to the varying wind speed and rotor speed. In order to select the most appropriate combination of Reynolds and Mach numbers for the various airfoils, the IEA 22 MW blade design was assumed as a reference, as it is the most recent and representative of the size and characteristics of upcoming offshore rotors. Individual airfoils are assigned to a representative span of the blade based on their thickness, as thicker airfoils are generally used in the inner part of the blade, while thinner ones are preferred in the outboard part to improve aerodynamic efficiency. Once the average chord of the blade section is computed, the inflow velocity for a given Reynolds number (Eq. 3) can be computed by inverting Eq. (1):

$$\begin{array}{}\text{(3)}& U={\displaystyle \frac{Re\mathit{\mu }}{\mathit{\rho }c}}.\end{array}$$

The chord length, Reynolds number, inflow velocity, and Mach number ranges assumed for each airfoil relative thickness are shown in Table 3. Five Reynolds numbers are simulated for each airfoil: 2.5×10⁶, 5×10⁶, 10×10⁶, 15×10⁶, and 20×10⁶. As a consequence of the smaller chord of the tip airfoils and since the same range of Reynolds number was simulated for all relative thicknesses, a maximum velocity of 153.3 m s⁻¹ is considered for airfoils with relative thickness of 21 % and below. This velocity should rather be regarded as a conservative upper bound and is unlikely to be reached during steady-state wind turbine operation. However, its inclusion ensures that the dataset remains valid under transient scenarios, such as rotor overspeed events, extreme gusts, or rapid changes in operating conditions, which can temporarily lead to significantly increased local inflow velocities. The portions of blade span upon which the airfoils are assumed to insist are shown in Fig. 3.

Table 3 Inflow conditions on the various airfoils based on relative thickness.

Figure 3 Assumed range of relative thickness as a function of blade span on the IEA 22 MW blade.

3 Methods

All data included in this dataset are computed using two-dimensional CFD. More specifically, compressible URANS simulations are performed in Ansys Fluent 23R2. The adoption of the URANS approach is not necessarily intended to capture unsteady flow features. Rather, time-dependent simulations are preferred over a steady-state approach as the former have shown to be less prone to numerical instabilities and more accurate, especially when laminar to turbulent boundary layer transition is modeled. The numerical set-up is based on a previously validated approach (Balduzzi et al., 2021; Papi et al., 2021; Venturi et al., 2024). The working fluid is dry air at sea-level pressure and 0 °C temperature (Table 2) and is modeled as an ideal gas.

An open-field, bullet-shaped domain is used as shown in Fig. 4. This geometry allows for the use of a single inlet boundary, where inflow velocity, turbulence intensity, and length scale are specified, and a single outlet boundary, where the free-stream pressure is imposed. The domain boundaries are placed 600 chords away to avoid interference with the pressure distribution on the airfoil, as discussed in more detail in Sect. 3.1, while also allowing the airfoil's wake to fully develop. The angle of attack is changed by rotating the airfoil inside the domain via a sliding interface. While not presently included in this dataset, this approach also allows for dynamic simulations to be performed. The k−ω SST turbulence closure (Menter, 1994) is used in fully turbulent simulations. The numerical decay of turbulence from the inlet to the airfoil is prevented using ad hoc source terms in the k and ω transport equations: this avoids artificially large length scales and consequent contamination of the solution in non-turbulent regions due to excessive eddy viscosity (Spalart and Rumsey, 2007). Boundary layer transition in the free transition simulations is modeled using the four-equation γ−Re_θ transition model (Langtry and Menter, 2009; Menter et al., 2004). In order to fine-tune the characteristics of the transition model, the correlation developed by Drela (1998) for the transition-critical Reynolds number is used, as explained in more detail in Sect. 3.2.

Figure 4 Domain used in computations.

The simulations are initialized from a steady Reynolds-averaged Navier–Stokes (RANS) simulation with the same inflow conditions and model set-up as the following URANS runs. RANS simulations are run for 2000 iterations or until all residuals drop below $\mathrm{1}\times {\mathrm{10}}^{-\mathrm{6}}$. After initialization, the simulations are run for 50 chord-based through-flow times or until numerical convergence, where numerical convergence is achieved if the percentage variation over the first and last values used to compute the mean aerodynamic coefficients is less than 0.5 %.

Experimental and numerical data generated during the AVATAR project (Schepers et al., 2018) serve as the benchmark for model validation and calibration in this study. During the project, experiments were conducted on the DU-00-W-212 airfoil at Reynolds numbers ranging from 3×10⁶ to 15×10⁶. Such high-Reynolds-number data are exceptionally rare and were made possible by the DNW High Pressure Wind Tunnel (HDG) in Göttingen. By pressurizing the tunnel to increase air density, flight-scale Reynolds numbers were achieved on a model with a chord of only 15 cm. The resulting dataset includes aerodynamic coefficients from positive to negative stall, alongside detailed surface pressure distributions. Furthermore, these benchmarks were used by project members to cross-verify various computational tools (Sørensen et al., 2016). The following sections detail the numerical framework: Sect. 3.1 discusses grid and domain sensitivity; Sect. 3.2 provides a verification of the fully turbulent set-up against established CFD results; Sect. 3.3 explores boundary layer transition modeling sensitivity; and Sect. 3.4 addresses the methodologies for post-stall extrapolation to high angles of attack.

3.1 Numerical model tuning and validation

Two-dimensional CFD computations of external aerodynamic bodies can be particularly sensitive to boundary conditions. To ensure truly open-field conditions, i.e., no artificial acceleration or distortion of the free stream, the boundaries need to be placed as far as 500 chords away from the airfoil (Golmirzaee and Wood, 2024). As highlighted by Sørensen et al. (2016), a smaller domain can be used if a point-vortex correction for inlet velocity is considered. While such a correction was not used in this study, it highlights the influence of lift on the flow-field distortion; thus large domains are required, particularly in high-lift situations. The influence of domain size on the predicted performance in fully turbulent conditions at a Reynolds number of 15×10⁶ is shown for the DU00-W-212 airfoil as a representative example in Fig. 5. Owing to the low Mach number of 0.08, the flow is treated as incompressible. The same grid resolution is maintained across all domain sizes by adjusting the total number of elements accordingly, thereby isolating the effect of domain size. A dedicated sensitivity analysis to the grid resolution is performed in the following step. As shown in Fig. 5, boundary distance particularly influences the estimation of the drag coefficient, and it is important to correctly capture the peak lift-to-drag ratio. In particular, as already suggested during the comparison rounds of the AVATAR project (Ozlem et al., 2017), the thresholds usually found in most of the literature studies (around 30 chords, Golmirzaee and Wood, 2024) are apparently not sufficient to ensure robust results. In the present study, the value of 600c (reported in solid lines in the graphs) has been adopted.

Figure 5 Relative error in maximum aerodynamic efficiency with respect to the largest tested domain (a), lift coefficient as a function of drag coefficient (b), and lift-to-drag ratio as a function of the angle of attack (c) with various domain sizes. Domain size is expressed as a function of the chord (C).

Focusing on the grid structure, a triangular unstructured grid combined with a prismatic boundary layer mesh is used in this study. The prismatic boundary layer mesh features 85 layers moving away from the airfoil's surface. This guarantees that the boundary layer is fully contained inside the prismatic layer in all tested conditions. The height of the first cell on the blade surface ensures a y⁺ value below 1 for all cases. The final mesh used in the validation of the numerical approach on the DU00-W-212 airfoil is shown in Fig. 6. Insensitivity to the mesh was ensured by comparing results obtained in free transition conditions at a Reynolds number of 15×10⁶ of the DU00-W-212 airfoil (Fig. 7). Free transition computations are considered more demanding from this point of view, as correctly capturing the transition point can require finer chordwise resolution than a fully turbulent computation (Menter et al., 2021). The free-stream mesh in the vicinity of the airfoil is parameterized depending on the number of elements along the airfoil itself; therefore increasing the number of elements along the airfoil also increases the total number of mesh elements in the domain, as shown in Table 4. To obtain grid-insensitive results, especially between 5 and 12° of angle of attack, the M3 grid was used in the subsequent simulations.

Figure 6 Triangular mesh used in computations. (a) Wake refinement and rotating region around the blade, (b) blade refinement, and (c) boundary layer mesh.

Table 4 Total and chordwise number of elements for tested computational grids.

Figure 7 Lift coefficient (a), drag coefficient (b), and lift-to-drag ratio (c) as a function of the angle of attack with varying mesh sizes. Reynolds number of 15×10⁶ and turbulence intensity of 0.34 %. Mesh size is expressed as a function of the number of chordwise elements (Ec). Reference data from Ozlem et al. (2017).

Figure 8 Lift coefficient (a), drag coefficient (b), and lift-to-drag ratio (c) as a function of the angle of attack with varying time step sizes. Reynolds number of 15×10⁶ and turbulence intensity of 0.34 %. Time step is expressed as a fraction of the through-flow time (Eq. 4). Reference data from Ozlem et al. (2017).

Sensitivity to the time step length is shown in Fig. 8. Although an implicit coupled solver is used to numerically solve the Navier–Stokes equations, which ensures numerical stability even for very high time steps in problems such as the one evaluated herein (Ansys Inc., 2024), the time step must be small enough to accurately resolve the relevant flow features. The time step is chosen as a function of the flow through-flow time, expressed as the ratio between the airfoil chord and the free-stream velocity as (Eq. 4)

$$\begin{array}{}\text{(4)}& tf=c/U.\end{array}$$

While lift coefficient is mostly insensitive to the chosen time step, the drag coefficient decreases if the time step is decreased from 5 to 10 time steps for chord through-flow time. A further decrease in time step does not lead to additional improvement. Based on this analysis, 20 time steps for chord through-flow times are used in the computation of the database.

3.2 Verification of the numerical model

The numerical set-up was verified with respect to the numerical simulations of the DU00-W-212 airfoil performed by Sørensen et al. (2016). The lift-to-drag ratio for the DU00W-212 airfoil at a Reynolds number of 3×10⁶ and a Mach number of 0.075 is shown in Fig. 9a. The pressure and skin friction coefficients at an angle of attack of 4°, obtained through an incompressible fully turbulent simulation of the DU00-W-212 in the same flow conditions, are shown in Fig. 9c and d. The data are overlaid on the numerical simulations from Sørensen et al. (2016), where the flow is considered incompressible in all the simulations. In addition to highlighting very good agreement between the current set-up and the reference data, Fig. 9a also shows a 2.2 % decrease in the lift-to-drag ratio at the angle of attack of peak efficiency of 8°. This measurable reduction in efficiency is present despite the relatively low inflow Mach number of 0.075, motivating the inclusion of compressibility in the dataset.

Figure 9 Lift-to-drag ratio as a function of the angle of attack (a) and pressure coefficient (b) and friction coefficient (c) as functions of the chordwise position. Reynolds number of 15×10⁶, turbulence intensity of 0.34 %, and Mach number of 0.075. Reference data in red from Sørensen et al. (2016) show the spread in model prediction that was found during the AVATAR project.

3.3 Boundary layer transition modeling

Laminar-to-turbulent transition is modeled with the γ−Re_θ formulation, which uses two additional transport equations, one for the intermittency γ and one for the momentum thickness Reynolds number Re_θ. The model relies on empirical correlations for the critical Reynolds number (Re_θ,t) and transition length, which are used to trigger turbulence production in the laminar boundary layer. Thus, they directly influence the position of the transition point and extension of the transition region. Despite being capable of modeling all kinds of boundary layer transitions, the original correlations (Langtry and Menter, 2009) were developed with turbomachinery flows in mind and were optimized for bypass transition occurring at free-stream turbulence values higher than 1 %. Despite working inside the ABL, modern wind turbine airfoils operate at notably lower turbulence levels, as the relatively steady tangential velocity due to the turbine's rotation is the main contributor to the local inflow, rather than the unsteady incoming wind. Such low-free-stream-turbulence inflows are dominated by the natural transition mechanism (Lobo et al., 2023). When simulating cases with low inflow turbulence, some authors have noted some deficiencies in the γ−Re_θ model (Diakakis et al., 2019). However, other studies have also shown that correlation with experiments can be significantly improved in correlation-based transition models through calibration (Colonia et al., 2016, 2017). For this reason, in this work the correlation developed by Drela (1998) for the critical Reynolds number Re_θt is used (Eq. 5):

$$\begin{array}{}\text{(5)}& \begin{array}{rl}R{e}_{\mathit{\theta }t}& \left(H,N_{\mathrm{cr}}\right)=\mathrm{155}\\ & +\mathrm{89}\left[\mathrm{0.25}\cdot \tanh \left({\displaystyle \frac{\mathrm{10}}{H-\mathrm{1}}}-\mathrm{5.5}\right)+\mathrm{1}\right]N_{\mathrm{cr}}^{\mathrm{1.25}}.\end{array}\end{array}$$

The correlation depends on the amplification factor N_cr and the boundary layer shape factor H. The former can be modeled as a function of the free-stream turbulence intensity (Eq. 6):

$$\begin{array}{}\text{(6)}& N_{\mathrm{cr}}=f\left(ti\right)=-\mathrm{8.43}-\mathrm{2.4}\ln \left({\displaystyle \frac{\mathrm{2.1}\tanh \left(\frac{\mathrm{ti}}{\mathrm{2.7}}\right)}{\mathrm{100}}}\right),\end{array}$$

or used as a calibration factor, as done in this study. The boundary layer shape factor, i.e., the ratio between boundary layer displacement thickness and momentum thickness, cannot be computed directly, as the γ−Re_θ transition model is based on local variables only and must be approximated with an empirical correlation. In this work, the correlation developed by White and Madani (2022), based on the experiments of Thwaites (1949) on flat plates under a streamwise pressure gradient, was used (Eq. 7):

$$\begin{array}{}\text{(7)}& \begin{array}{rl}H\left({\mathit{\lambda }}_{\mathit{\theta }}\right)& =\mathrm{2}+\mathrm{4.14}x-\mathrm{83.5}{x}^{\mathrm{2}}+\mathrm{854}{x}^{\mathrm{3}}-\mathrm{3337}{x}^{\mathrm{4}}+\mathrm{4576}{x}^{\mathrm{5}}\\ & x=\mathrm{0.25}-{\mathit{\lambda }}_{\mathit{\theta }},\end{array}\end{array}$$

where λ_θ is the local pressure gradient, computed locally from the resolved flow field. Once again, the approach was validated with respect to high-Reynolds experimental tests from the AVATAR project (Pires et al., 2016). Results obtained using the γ−Re_θ transition model, with the correlation for Re_θt(H,N_cr) described in Eq. (4), are compared to results obtained with the same transition model with the default correlation developed by Langtry and Menter (2009) and to experimental data in Fig. 10. Although the dataset contains measurements for various Reynolds numbers, ranging from 3×10⁶ to 15×10⁶, the tests differ not only in Reynolds number but also in turbulence intensity. Therefore, the test at Re 15×10⁶ was chosen for calibration of the numerical models as this condition is the most representative of very large offshore blades. As shown in Fig. 10, while the default correlation was found to produce good results in the conditions tested herein, using the correlation developed by Drela (Eq. 5) with N_cr=6 reduced the discrepancy in the lift-to-drag ratio at low angles of attack and was thus chosen for the computation of the entire dataset. The free-stream turbulence intensity of the simulations is set to 0.34 %, which matches the AVATAR experiments at Re15×10⁶.

Figure 10 Lift coefficient (a), drag coefficient (b), and lift-to-drag ratio (c) as a function of the angle of attack. Experiments from the AVATAR project (Pires et al., 2016) are compared to the results obtained with the numerical set-up described in Sect. 2 and various correlations for the γ−Re_θ transition model on the DU00-W-212 airfoil. Reynolds number of 15×10⁶, turbulence intensity of 0.34 %, and reference data from Ozlem et al. (2017).

3.4 Post-stall extrapolation

As per common practice, the raw lift and drag coefficients were first manually cleaned by removing inaccurate/unconverged results, generally close to or after the positive or negative stall points. Clean coefficients were then extrapolated to the full ±180° angle of attack range using empirical post-stall correlations. Differently from other studies, special attention was given to the extrapolation process. In particular, the choice and calibration of the post-stall extrapolation method were driven by dedicated scale-resolving stress-blended eddy simulations (SBESs), which were carried out for the FFA-W3-241 airfoil at a Reynolds number of 10×10⁶ and selected angles of attack (see Fig. 11).

Figure 11 Q-criterion isosurfaces for the FFA-W3-241 airfoil at 30° (a, c) and 90° (b, d) angle of attack from scale-resolving SBESs. Reynolds number of 10×10⁶. Side view (a, b) and top view (c, d).

SBES blends an LES turbulence closure in the free stream and wake of the airfoil when the grid resolution allows it and a URANS approach in the boundary layer region (Menter, 2021). Simulations are performed in the Ansys Fluent 24R1 GPU-accelerated solver using the k−ω SST turbulence closure in the airfoil boundary layer and the WALE LES sub-grid-scale model in the free stream. As shown in Fig. 11, the airfoil section is extruded six chords to allow for spanwise flow non-uniformities to fully develop and avoid the “lock-in” of wake vortex structures and the consequent overestimation of force fluctuations (Mittal and Balachandar, 1995). The same bullet-shaped domain as used for the two-dimensional computations is also used in the three-dimensional case, although the inlet and outlet are placed 70 chords away from the airfoil to avoid low-quality cells due to an excessive aspect ratio. Periodic boundary conditions are imposed on the sidewalls of the three-dimensional section, thereby ensuring an infinite aspect ratio in the simulation and the absence of tip effects.

This set-up is compliant with the recommendations found in CFD prediction of airfoil deep-stall performance using improved delayed detached-eddy simulations Sørensen and Timmer (2017), where at least four chords of spanwise distance are recommended for simulations with a high angle of attack. The boundary layer is modeled with a series of prismatic elements, ensuring y⁺ values at the airfoil surface below 1. The mesh in the immediate vicinity of the airfoil is also prismatic, while a polyhedral mesh is used in the airfoil wake and is chosen due to generally lower numerical diffusion and improved gradient reproduction with respect to tetrahedral cells (Wang et al., 2021). A total of 400 elements are used along the airfoil in both the chordwise and the spanwise directions. This choice is consistent with the state-of-the-art of scale-resolving simulations (Manolesos and Papadakis, 2021; Sørensen and Timmer, 2017) in the chordwise direction, while the total number of elements in the spanwise direction is limited due to solver memory constraints. Due to the expected occurrence of large-scale, strongly three-dimensional flow structures in the deep-stall cases analyzed herein, priority is given to the total resolved blade span rather than to the spanwise discretization. Nevertheless, an aspect ratio of approximately 1.5–1.8 is achieved for the cells in the airfoil's proximity, quickly transitioning to lower values in the polyhedral wake region. The total number of elements is roughly 24×10⁶. Similarly to CFD prediction of airfoil deep-stall performance using improved delayed detached-eddy simulations (Sørensen and Timmer, 2017), the domain extends 50 chords from the airfoil and is bullet-shaped (see Fig. 4). Simulations are first initialized with a RANS simulation and then run for 300 chord-based through-flow times (Eq. 8):

$$\begin{array}{}\text{(8)}& \mathit{\tau }=t\cdot U/c,\end{array}$$

or until numerical convergence of the mean lift and drag coefficients. Results are considered converged when the mean values of the aerodynamic coefficients vary by less than 2 % over the final 70 τ. Lift and drag coefficients are then averaged starting from τ= 60. The SIMPLE pressure–velocity coupling solution algorithm is used. The time step is chosen in order to guarantee a Courant number (Courant et al., 1967) that does not exceed 2 in most cells around the airfoil and is below 5 in the entire domain. This guarantees the stability and robust convergence of the SIMPLE algorithm. The strongly three-dimensional characteristics of the flow at both 30 and 90° of flow incidence in particular can be seen in Fig. 11. They are most noticeable at 90° angle of attack (Fig. 11b, d), where a large-scale vortex structure can be seen in the airfoil's wake. The total computational cost for SBESs is approximately 1400 GPU hours using a single NVidia A100 80 GB GPU.

After initial testing of various post-stall extrapolation models, such as those developed by Viterna and Janetzke (1982), Spera (2008), Beans and Jakubowski (1983), and Kirke (1998), the choice was narrowed down to the models developed by Montgomerie (2004) and more recently by Battisti et al. (2020), based on their agreement with scale-resolving simulations (Fig. 12). Based on the comparison of the two methods with SBESs of the FFA-W3-241 airfoil (Fig. 12), drag is extrapolated to ±180° using the method by Battisti et al. (2020), while the Montgomerie (2004) method is used for lift. While a direct comparison to high-angle-of-attack experimental data was not directly performed in this study, the maximum drag predicted by the SBES is in line with the existing literature, towards which the reader is redirected for more in-depth discussions on the topic (Timmer, 2010, 2020). The moment coefficient is extrapolated using the relation in Laino and Hansen (2002).

Figure 12 Lift (a) and drag (b) coefficients for the FFA-W3-241 airfoil at a Reynolds number of 10×10⁶. Two-dimensional CFD simulations extrapolated with the methods proposed by Montgomerie (2004) and Battisti et al. (2020) compared to scale-resolving SBESs.

After post-stall extrapolation, fully turbulent and free transition coefficients are blended with a 70 % free transition and 30 % fully turbulent ratio. These sets of coefficients can be considered representative of a new, clean blade, while also including the effect of aerodynamic deterioration which may occur during operation. This practice is fairly common in academia and has been followed in the design of the latest-generation offshore reference rotor, the IEA 22 MW RWT (Zahle et al., 2024).

3.5 Static parameters

In a final step, the static parameters used in the Beddoes–Leishman (BL)-type dynamic stall models underlying most aero-elastic tools (Leishman, 2016) are computed from the obtained polar data. A detailed description of each parameter can be found in Table 5. Although not strictly related to the airfoil dynamic behavior, these coefficients represent two-thirds of the parameters necessary to calibrate the model and have a strong influence on the overall accuracy (Melani et al., 2024). The resulting coefficients are provided in the database in AeroDyn format (Jonkman et al., 2025).

Table 5 Static parameters for Beddoes–Leishman (BL)-type dynamic stall models used in the current study.

The parameter extraction follows the workflow outlined in Melani et al. (2024). For the sake of brevity, the process is described only for positive angles of attack, but the same considerations apply to the negative ones. As a preliminary step, lift (C_l) and drag (C_d) coefficients are converted to the normal (C_n) and chordwise (C_c) ones used in the BL formulation. Parameters governing the airfoil behavior in the attached- and separated-flow regimes are then determined sequentially.

Figure 13 Visual representation of the computation of the different parameters required by AeroDyn dynamic stall models from static polar data: (a–c) attached-flow coefficients and (d–f) separated-flow coefficients.

For the attached-flow regime, constants α₀, C_Nα, C_D0, C_M0, and η are directly obtained from the C_n, C_c, and C_m curves. The zero-lift angle of attack α₀ is located at the intersection of the C_n(α) curve with the α axis, while C_M0 and C_D0 are the moment and drag values at α=α₀ (see Fig. 13a–c). The curve slope C_Nα is evaluated via linear regression of the data in the linear region, delimited in Fig. 13a by the interval [${\mathit{\alpha }}_{\mathrm{lin},\mathrm{neg}}{\mathit{\alpha }}_{\mathrm{lin},\mathrm{pos}}$]. This accounts for possible noise in the static data. In this work, an innovative automated method to identify this region is adopted: (i) progressively narrower intervals, centered on α₀, are considered; (ii) for each interval, data are fitted via linear regression, using the corresponding R² as a linearity index; and (iii) the linear region is identified as the region where R² > 0.99 at both bounds. Finally, the recovery factor η_e is selected by minimizing the discrepancy between ${\mathit{\eta }}_e{C}_{\mathrm{c}}^{\mathrm{pot}}$ (see Table 5) and the pressure component of chordwise force $C_{\mathrm{c}}^{\mathrm{stat}}=C_{\mathrm{l}}\sin \mathit{\alpha }-(C_{\mathrm{d}}-C_{\mathrm{D},\mathrm{0}})\cos \mathit{\alpha }$ in the linear region (see Fig. 13b).

In the BL model, flow separation is described by the static separation point (f) and center of pressure (x_CP). The separation point is derived from the static characteristics via the Kirchhoff formula in Eq. (9):

$$\begin{array}{}\text{(9)}& f\left(\mathit{\alpha }\right)={\left(\mathrm{2}\sqrt{{\displaystyle \frac{C_{\mathrm{n}}}{C_{N\mathit{\alpha }}\left(\mathit{\alpha }-{\mathit{\alpha }}_{\mathrm{0}}\right)}}}-\mathrm{1}\right)}^{\mathrm{2}}.\end{array}$$

The resulting f(α) curve is fitted with a piecewise exponential (Leishman) function, reported in Table 5, yielding the parameters α₁ (f= 0.7), S₁, and S₂. The process is illustrated in Fig. 13d. The center of pressure is obtained in a similar fashion from the knowledge of C_n and C_m, as in Eq. (10):

$$\begin{array}{}\text{(10)}& x_{\mathrm{CP}}\left(\mathit{\alpha }\right)={\displaystyle \frac{C_m-C_{M\mathrm{0}}}{C_{\mathrm{n}}}},\end{array}$$

and subsequently fitted with the semi-empirical formulation of Beddoes (see Table 5). The procedure, as well as the meaning of the different fitting constants k₀, k₁, k₂, and k₃, is shown in Fig. 13f.

Finally, the static stall threshold C_N1, used in the BL model to mark leading-edge separation, is identified as the C_n value at the maximum C_c in the pre-stall region. This criterion is deemed to be the most robust among the criteria available in the literature, as the C_c curve always presents one peak at the static stall point. The latter is also more coherent with the algorithm used in the BL model to reconstruct the C_c curve (see Table 5), where C_N1 is the breakpoint between the formulations used for the pre- and post-stall regions (see Fig. 13e). The shape of the post-stall curve can be further tuned by adjusting the parameters D_CC and ${\widehat{k}}_{\mathrm{1}}$, though D_CC is not used in AeroDyn.

3.6 Application of database to wind turbine simulations

As discussed in Sect. 2, this dataset is intended for use in wind turbine applications. Therefore, it must be noted that the coefficients in this dataset are relative to a non-rotating infinite-span uniform wing. In other words, they are relative to the two-dimensional sections described in Sect. 2. Before being applied to a three-dimensional wind turbine blade, they should be corrected for three-dimensional rotational augmentation (Snel et al., 1992; Timmer, 2010) and for finite aspect ratio, which affects the post-stall extrapolation and maximum drag coefficient (Battisti et al., 2020; Timmer, 2010; Viterna and Corrigan, 1982).

4 Annotated comparisons of airfoil families

The operating condition of an airfoil, often expressed in terms of Reynolds and Mach numbers, greatly influences its performance. While this is commonly acknowledged and is not a novelty in itself, the variety of conditions of the examined dataset, which have been tailored specifically to multi-megawatt wind turbine rotors, nevertheless enables the identification of trends directly relevant to rotor aerodynamic design and optimization.

As a first example, lift and drag coefficients, as well as the aerodynamic efficiency, are shown as a function of the angle of attack for the FFA-W3-241 airfoil in Fig. 14. As expected, the Reynolds number has a strong influence on all quantities of interest. For both the lift and the drag coefficients, fully turbulent simulations are more sensitive than transitional ones to the changes in inflow conditions, as a general decrease in drag coefficient with increasing Reynolds number can be noted. The stall limit (static stall angle and maximum lift coefficient before stall) increases instead. These changes are reflected in airfoil efficiency, which grows with the Reynolds number in fully turbulent conditions (Fig. 14c).

Figure 14 (a) Lift, (b) drag, and (c) aerodynamic efficiency of the FFA-W3-241 airfoil for various Reynolds numbers in fully turbulent and free transition conditions.

More general conclusions can be drawn from Figs. 15 and 16, where the design angle of attack, here taken as the point where the airfoil presents the maximum aerodynamic efficiency $E=C_{\mathrm{l}}/C_{\mathrm{d}}$, is shown together with the corresponding aerodynamic efficiency and lift coefficient as a function of Reynolds and Mach number. The same trends are reported for the static stall angle. These are considered the key parameters for assessing the performance of an aerodynamic section used in a wind turbine rotor. Lift-to-drag ratio is of paramount importance to minimize the penalty on torque (drag) per unit of useful aerodynamic force (lift). In a horizontal-axis wind turbine, this is especially important in the outboard sections of the blade, where, due to the blade's high tangential velocity, lift is directed mostly out of plane, thereby increasing unwanted thrust force rather than the driving torque, and drag is directed mostly in-plane, accentuating its detrimental effect on torque. The maximum airfoil efficiency for airfoils with 30 % or lower relative thickness in fully turbulent conditions is shown in Fig. 15 as a function of Reynolds number (c) and Mach number (g), respectively. For all airfoil families, lower thickness ratios generally lead to higher aerodynamic efficiency. For all airfoils, efficiency tends to increase with the Reynolds number. However, particularly for more outboard, thinner airfoils, efficiency plateaus with very little, if any, increase beyond 15×10⁶ Reynolds number. From an inspection of Fig. 15g, compressibility appears to influence this effect. In fact, while the Reynolds number is the same for all airfoils and ranges from 2.5×10⁶ to 2×10⁷, the Mach number is not, as thinner airfoils are simulated with higher inflow velocities. If maximum efficiency is shown as a function of the Mach number (Fig. 15g), with the only exception of the DU91-W2-250 and FX77-258 airfoils, which underperform compared to the rest, all airfoils follow a very similar trend. When compressibility effects are low and the efficiency curves remain below Ma=0.2 (as is the case for the FFA-W3-301 airfoil for example), efficiency increases as the Reynolds (and Mach) number increase. Above Ma=0.25 however, the maximum lift-to-drag ratio starts to plateau for all the tested airfoils, with the FFA-W3-211 airfoil being a clear example of this trend. Both the Reynolds number and compressibility have a marked effect on stall performance. In fact, both the stall angle (Fig. 15e) and the design angle (Fig. 15f) tend to initially increase at lower Reynolds numbers, due to the increased resistance to stall as the airfoil boundary layer becomes increasingly thinner, but tend to decrease at higher Mach numbers, highlighting the effect of compressibility.

Figure 15 Stall angle computed as angle of attack of maximum lift (a, e), design angle of attack (b, f), maximum aerodynamic efficiency (c, g), and lift at design angle of attack (d, h). Quantities shown as a function of the Reynolds number (a–d) and Mach number (e–h). Fully turbulent simulations of airfoils with less than 30 % relative thickness.

Figure 16 Stall angle computed as angle of attack of maximum lift (a, e), design angle of attack (b, f), maximum aerodynamic efficiency (c, g), and lift at design angle of attack (d, h). Quantities shown as a function of the Reynolds number (a–d) and Mach number (e–h). Free transition simulations of airfoils with less than 30 % relative thickness.

The variation in the design lift coefficient as a function of the Reynolds and Mach number is also significant, with maximum differences exceeding 20 % in some cases. The differences are stronger for more inboard airfoils, which are less affected by flow compressibility. These variations may have a significant impact on design, as blade sections can experience significant variations in Reynolds and Mach number depending on the rotational speed.

Figure 17 (a) Pressure coefficient and (b) skin friction for the FFA-W3-241 airfoil in fully turbulent conditions (blue) and free transition conditions (red) at an angle of attack of 8°.

The same analysis is repeated in free transition conditions. The stall angle, the design angle of attack, the lift coefficient at the design angle of attack, and the maximum airfoil efficiency are shown as a function of the Reynolds and Mach numbers for various airfoil families in Fig. 16. While the effect of compressibility on the stall angle, which generally tends to be anticipated with increasing Mach number, is similar to what was noted in fully turbulent conditions, the variations in design angle of attack, maximum efficiency, and design lift coefficient follow different trends. When the boundary layer is allowed to transition, the model predicts greater resistance to stall, as demonstrated by the higher stall angles, especially at lower angles of attack. Resistance to stall is less influenced by the Reynolds number and appears to be more airfoil dependent. In addition, aerodynamic efficiency follows a different trend as a function of Reynolds number with respect to what was noted in fully turbulent conditions (Fig. 15c) and generally decreases as Reynolds number increases. This trend is also confirmed by experiments, although inflow characteristics, and particularly inflow turbulence, varied slightly in the AVATAR test campaign (Pires et al., 2016). The reason for the decrease in maximum efficiency can be explained by examining Fig. 16, where the pressure (C_p) and skin friction (C_f) coefficients for the FFA-W3-241 airfoil in both free transition and fully turbulent conditions are shown. In both cases, pressure loading in the front region of the airfoil suction side ($x/c\mathit{<}\mathrm{0.4}$) tends to increase with Reynolds number due to a reduction in boundary layer thickness. In fully turbulent conditions, this increase is met by a corresponding decrease in skin friction, explaining the increased maximum efficiency. When the boundary layer is allowed to transition from laminar to turbulent however, the higher the Reynolds number, the earlier the flow tends to transition, as shown by the increases in skin friction in Fig. 17b. The earlier laminar-to-turbulent transition reduces the low-friction area of the airfoil operating with a laminar boundary layer. The anticipated transition point is balanced by the decrease skin friction in both the laminar and the turbulent regimes at higher Reynolds numbers (Fig. 17b), ultimately leading to the trends shown in Fig. 16c.

Figure 18 Comparison of simulated trends (gray lines) in free transition conditions with experimental results of the DU-00-W-212 airfoil from Ozlem et al. (2017) (black markers). Design angle of attack (a) and maximum aerodynamic efficiency (b).

A similar trend in terms of the variation in maximum efficiency and design angle of attack as a function of the Reynolds number is also noted in the experiments performed by Ozlem et al. (2017) on the DU-00-W-212 airfoil. The experimental results are compared in Fig. 18 to the same set of airfoils shown in Fig. 16. In the experiments performed on the DU-00-W-212 airfoils, efficiency peaks at $Re=\mathrm{9}\times {\mathrm{10}}^{\mathrm{6}}$ (Fig. 18b), while the angle of attack at which the maximum lift-to-drag ratio is located (Fig. 18a) decreases as the Reynolds number increases before starting to increase at $Re=\mathrm{1.5}\times {\mathrm{10}}^{\mathrm{7}}$. The simulations are consistent with this trend, with some differences from airfoil to airfoil and depending on the inflow conditions. The Mach number is also much higher in the simulations, as it ranges between 0.01 and 0.46, while it does not exceed 0.1 in the experimental reference.

The absolute differences between the fully turbulent and free transition cases in terms of stall and design angle of attack, aerodynamic efficiency, and design lift coefficient are shown in Fig. 19. Small differences between the performance in fully turbulent and in free transition conditions are generally desirable. These two extreme conditions lie at either end of the operational windows of a wind turbine blade and correspond to a dirty, worn-down blade or a new, clean blade, respectively. The larger the difference in terms of efficiency between these two configurations, the larger the performance degradation during the turbine lifetime will be. As shown in Fig. 19, there are generally larger differences between the free transition and fully turbulent scenarios at lower Reynolds numbers. This trend is most noticeable when examining maximum efficiency, as the gap between the free transition and fully turbulent cases closes as the fully turbulent efficiency increases.

Figure 19 Difference between values of free transition and fully turbulent simulations. Stall angle computed as angle of attack of maximum lift (a, e), design angle of attack (b, f), maximum aerodynamic efficiency (c, g), and lift at design angle of attack (d, h). Quantities shown as a function of the Reynolds number (a–d) and Mach number (e–h).

Figure 20 Angles of attack where efficiency exceeds 85 % of the maximum value (lines), angle of attack of maximum efficiency (filled markers), and stall angle of attack (cross) for various airfoils. The ranges for the angle of attack are shown for the five Reynolds numbers in the dataset. The y values correspond to the maximum efficiency. 21 % thick airfoils (a, d), 24 %–25 % thick airfoils (b, e), and 30 % thick airfoils (c, f). Fully turbulent (a, b, c) and free transition (d, e, f) conditions.

The operating range of some of the airfoils in the database is compared in Fig. 20 by showing the ranges of angles of attack where efficiency exceeds 85 % of the maximum value at each Reynolds number. A broad high-efficiency operating range is desirable in wind turbine applications as these machines operate in highly turbulent inflow and face highly variable inflow conditions, which cause the angle of attack to vary during operation. In addition, modern wind turbine blades tend to deform as rotor load increases due to combined effects of aerodynamic loading and aeroelastic phenomena such as bend–twist and shear–twist coupling. As such, a broader operating range can increase design options and the optimal operational window. When fully turbulent and free transition conditions are compared (top and bottom rows of Fig. 20, respectively), the range of angles of attack with good performance remains similar, but the stall margin, i.e., the distance of this optimal operating range from the static stall point, is higher in the case of a clean blade. In free transition conditions, however, some airfoils, such as the DU08-W-210, OSO-21, and DU91-W2-250, suffer from reductions in the optimal operating window at high Reynolds number due to the effects of compressibility, which leads to early stall and performance degradation. FFA-W3 airfoils generally show good performance in this regard.

5 Conclusions

This paper describes a database of wind turbine airfoil aerodynamic coefficients. The aerodynamic coefficients are computed using a state-of-the-art CFD methodology, which has been fine-tuned and validated based on available experimental data. Results include the effects of compressibility, and computations are performed with and without laminar-to-turbulent transition. The total computational resources required to compute the database exceed 100 000 CPU hours on modern hardware (dual AMD-EPYC 7413 server). This dataset expands on measurement-based resources such as the Aerodynamic Table Generator developed by TNO, which is constrained by limited Reynolds number coverage, and provides a high-fidelity complement to similar recent datasets (Karcher et al., 2025) generated via viscous panel methods. The quantity and fidelity of the generated data allow for the effects of Reynolds and Mach number variation to be investigated across different families of wind turbine airfoils featuring relative thicknesses between 21 % and 36 %. In particular, these effects are more pronounced when considering fully turbulent computations rather than in the free transition cases. In fact, when not considering boundary layer transition, peak efficiency increases as the Reynolds number increases, but benefits are greater at lower Mach numbers. Above a free-stream Mach number of 0.2, efficiency gains appear to be limited by the effects of compressibility. In free transition conditions, the variation in the transition point due to the change in inflow conditions can compensate for the effects of Reynolds and Mach number variation, and the variation in airfoil efficiency as a function of these non-dimensional parameters is much more case dependent. Because of these trends, the difference in efficiency between the free transition case, representative of a clean blade, and the fully turbulent case, representative of a soiled blade, decreases as the Reynolds number increases. Based on this evidence, the benefits of operating a clean blade are more strongly felt at lower wind speeds, where the Reynolds number is lower, rather than at higher wind speeds, where the Reynolds number increases due to the higher rotational speed. Finally, the more recent airfoils included in this database, namely the OSO family, outperform other academic open-access families, such as the FFA-W3 and DU airfoils, over a large range of relative thicknesses, proving that modern design tools can indeed improve blade sectional performance. On the other hand, the variations in design angle of attack as a function of the Reynolds number and boundary layer modeling (free transition compared to fully turbulent) are more significant than the OSO family of airfoils, especially for the geometries with higher relative thickness. Finally, it is worth remarking that the presented dataset is to be intended as a “living” dataset, and new geometries may be added to it in the future based on the wind energy community's requests.