Sensitivity of Dynamic Stall Models to Dynamic Excitation on Large Flexible Wind Turbine Blades in Edgewise Vibrations

Abstract

Present studies are specifically aimed at investigating the sensitivity of different dynamic stall models when exposed to various excitation frequencies. The investigations are targeted at blade edgewise vibrations. This is carried out on a modified version of the IEA 15 MW reference wind turbine employing a wind turbine design tool, DNV Bladed. State-of-the-art dynamic stall models for wind turbine applications, such as the Øye model, Beddoes–Leishman (BL) model and the newly developed IAG model, are evaluated. The beginning of the research work starts by evaluating different dynamic stall models on rigid blade section forces against known airfoil datasets. Then, the blade flexibility is considered to enable systematic evaluations of the blade flexibility influences in comparison to the rigid blade cases. It is observed that the range of the angle of attack grows depending on the excitation frequency and the adopted dynamic stall model. The critical excitation frequency range and the effects of twist distribution are then identified from the studies, which can be useful as a rough guidance when designing wind turbine blades.

1. Introduction

Wind turbine loads are heavily influenced by the operating conditions and environments [1]. Therefore, standards [2,3,4,5] state that blade designers/manufacturers need to run thousands of load cases that simulate wind turbines under various scenarios in the design and certification processes. The test cases range from normal power production (typically aimed at fatigue load analysis) to extreme operations when the turbine is shut down, either in parked or idling modes. The extreme cases are typically applied for high wind speeds above the cut-out wind speed of the turbine. In such a circumstance, the wind turbine controller is normally no longer active when no backup power supply is available. Hence, the blade is prone to operating at large angle of attack depending on the yaw angle, and it remains locked in this position. This could lead to strong unsteady vibrations and the aeroelastic behavior of the blade can be dangerous to its own structural integrity.

Despite the clear environmental challenges, wind turbine blades must survive even in harsh a climate, as defined in the set of loads calculation within the IEC [2] and DNV [3] standards. The consideration of unsteady aerodynamic effects has become even more important in the last few decades since the blade size has increased considerably, making the blade more slender and flexible. The aeroelastic characteristics of such long blades are not well understood. Turbine instability and vibrations often become huge concerns and limiting factors for designs.

The wind turbine aeroelastic effects are strongly affected by the amount of damping on the blade itself, which can generally be divided into structural damping and aerodynamic damping. Structural damping represents the structural design and the material used to manufacture the blade itself. Frequency response analysis is commonly used to estimate the structural damping on composite materials [6,7,8]. However, predicting the aerodynamic damping is not a trivial task because it can significantly change depending on the operating conditions of the turbine. This type of damping is affected by the acting forces on the blade, which alter the angle of attack due to elastic response. When the total damping of the blade becomes negative, instability of the blade is enhanced. It is well understood that the edgewise direction is prone to negative damping.

The aeroelastic phenomenon in the extreme situation is highly unsteady, and as a direct consequence, the induction-based blade element momentum theory commonly adopted in engineering tools shall be corrected using a dynamic stall model. Examples of such models are the Beddoes–Leishman (BL) model [9], Larsen model [10], Øye model [11], Snel model [12], Risø model [13], IAG model [14,15], ONERA model [16], Boeing–Vertol model [17], etc. Having a clear picture of how these dynamic stall models behave in challenging operations will be an important factor for quantifying the impacts on aeroelastic characteristics.

Several studies have attempted to evaluate the dynamic stall characteristics of wind turbines in standstill conditions (parked/idling) [18,19,20,21,22,23,24,25], and it is commonly believed that dynamic stall models tend to underestimate the true aerodynamic damping due to inaccurate hysteresis computations [26]. Hansen [19] provided a proposal for an improved modal analysis to improve the understanding of wind turbine stability. Sarkar and Bijl [27] evaluated non-linear aeroelastic behavior of an oscillating airfoil during stall-induced vibration in two-dimensional environments. Wang et al. [20], for instance, focused their studies on the linear model predictions when using the BL model and the ONERA model. Bir and Jonkman [28] investigated the aeroelastic instability of a floating 5 MW wind turbine, and pointed out the distinction between aeroelastic instability and resonance, where the latter requires external excitation at a distinct frequency. Tibaldi et al. [29] performed a numerical investigation on wind turbine resonant vibrations, where it was highlighted that depending on the external excitation, different aeroelastic modes can be excited. The responses of different dynamic stall models for a wind turbine in standstill conditions are further evaluated, for example, by Wang et al. [20], Faber [30], Bangga et al. [18], Andreou [31], Bangga and Yu [22]. In Ref. [32], commonly adopted engineering approaches were reviewed and discussed. Horcas et al. [33] provided an improved understanding of standstill flow by adopting high-fidelity computational fluid dynamics (CFD) approaches, but direct adoption of the studies to load case analysis is still impractical due to the excessive computational effort.

The present paper is specifically aimed at addressing the critical issue of how sensitive dynamic stall models actually are when they are exposed to dynamic excitation in the edgewise direction. To the best of the author’s knowledge, comprehensive studies in this area are lacking for extremely large wind turbine blade structures. This is especially crucial when different dynamic stall models and excitation frequencies are considered, leading to uncertainties in wind turbine blade design and certification processes. Existing references [18,19,20,21,22,23,24,25,27,29,30,31,32,32] were either limited to small to medium wind turbine size or were focused on generic responses of the wind turbine without detailed examinations on the sectional characteristics. No studies have been performed to specifically analyze the dynamic stall model responses to dynamic excitation in edgewise vibrations. Considering the current state-of-the-art and the research gap, the present paper investigates the aeroelastic sensitivity of a large wind turbine blade using three different dynamic stall models commonly applied in the designs, namely (1) the Øye model [11], (2) the incompressible version of the BL model [3] (Risø model [13]) and the newly developed IAG model [15]. The studies start by comparing different dynamic stall models against a pitching airfoil case, where experimental data are available. Then, the investigations are extended by exciting the rotor blade by means of an oscillating wind at various frequencies, ranging from below to above the first edgewise modal frequency of the blade. The growth of the angle of attack response is monitored at different wind speeds to quantify the direct influence of the unsteady aerodynamic models on the aeroelastic responses.

The present paper is organized as follows. Section 2 provides a description of the employed code, approaches and the test cases of interest. The results are discussed in Section 3.1, Section 3.2, Section 3.3 and Section 3.4, starting with the comparison of pitching airfoil data, then it continues with the aeroelastic response at the blade level. Finally, the paper will be concluded with suggestions for follow-up studies in Section 4.

2. Methodology

2.1. Wind Turbine Design Code Bladed

Bladed is an engineering tool used in numerous wind turbine industries for designs and certifications as well as in a research environment. The code solves aero-servo-hydro-elastic problems of wind turbines by utilizing blade element momentum (BEM) theory coupled with structural and control algorithms. The rotor aerodynamic solutions are obtained using time-marching approaches of the induced velocities at each blade section. Bladed uses a local BEM implementation [34]. The solutions are coupled with structural dynamic models, which equip a multibody dynamics approach. A multi-part blade method is adopted to model larger deflections using linear beam theory. Efficient solutions of these complex nonlinear problems are made possible via the fixed-step implicit Newmark-$\beta$ integrator. Several publications documented the verification studies of the Bladed code performance and accuracy [15,35,36,37,38,39]. Bladed version 4.16.0 is utilized in the present studies.

In the present studies, three dynamic stall models in Bladed are adopted, namely Øye [11], Beddoes–Leishman (BL) [34] and a recently developed IAG model [15]. The incompressible version of the BL model [34] has been widely adopted in wind turbine designs and can be considered as the state-of-the-art model. Although the model works sufficiently well for most conditions, the accuracy of the model is constantly challenged when it operates under complex inflow conditions with a large yaw misalignment. This implies that the angle of attack of the blade sections can be significantly large [15]. The IAG model has recently been developed specifically for wind turbine applications and has been tested against experimental data in various scenarios [14,15,40]. The robustness and applicability of the IAG dynamic stall model for wind energy field were further improved in the second-generation model by enabling a formulation in the state-space representation, as well as by removing the compressibility effects and by decreasing the number of required constants [15]. In Ref. [15], it was demonstrated that the IAG model performs well against experimental data and it outperforms other models when the airfoil operates at large angles of attack beyond stall, and the tests were carried out at various reduced frequencies and amplitudes. The IAG model was further compared against other dynamic stall models for various load cases in Refs. [18,22].

2.2. Blade Specification and Airfoil Data Treatment

The present studies mainly adopt a blade constructed by the S809 airfoil. This airfoil is chosen because there are experimental data available for both static and dynamic conditions [41]. Thus, validation of the results for the pitching airfoil case will be possible. This airfoil has a relative thickness of $t/c$ = 21% and was measured at a Reynolds number of a million. For the present studies, results based on the airfoil equipped with a leading edge grit (turbulator) will be adopted. This is performed to simulate the "soiled" effects on a wind turbine blade in the field.

In the present studies, the blade geometry of the IEA 15 MW reference wind turbine version 1.1 is adopted [42]. The blade is selected since it is representative of the current blade design that is associated with all complexity and highly flexible structures. The blade has a length of about 117 m. Starting from 91.33 m of the blade length, the relative thickness of the airfoil section (relative to the chord) is about 21%. In the present studies, the IEA 15 MW reference wind turbine blade is modified by replacing the applied airfoils at all blade sections. First, the airfoil data from 91.33 m of the blade length are replaced with the S809 airfoil. The remaining parts of the blade up to the root area are then generated based on interpolation (according to relative thickness). The root of the blade adopts a cylinder with a constant drag coefficient of about 0.35. This process is illustrated in Figure 1 and Figure 2. It is noted that most investigations will be focused on the dynamic response of the blade section at 93.68 m, except when stated otherwise. A more realistic test case can be designed by using the actual blade airfoils without modification. However, this will involve complex behavior of different airfoils with different stall characteristics. This makes it difficult to draw conclusions and to analyze further. Thus, the simplification performed in the present paper is used as a starting point. This simplification is deemed necessary to eliminate different stall characteristics of different airfoils, thus enabling the studies to focus on the dynamic response of the simplified blade. This way, the actual sensitivity of the dynamic stall model can be explored in isolation.

Since the S809 airfoil was measured only to a limited range of the angle of attack (α), it is a common practice to extrapolate the data to cover the α = −180° to α = +180° regime. In this paper, a modified Viterna approach [43] is adopted and the results are plotted in Figure 2. The improved method developed in the present studies gives an enhanced extrapolation accuracy for the pitching moment coefficient. A full description of the method is given in Appendix A. A 360° extrapolation is usually adopted in the early stages of wind turbine designs to give a rough idea about the turbine loads before any further detailed assessments are made. The usage of the extrapolated data will affect the accuracy of wind turbine calculations. The reconstruction accuracy deteriorates the shorter the data one has. It is highlighted in Ref. [36] how the dynamic calculations are influenced by the steady data reconstruction, and adopting inaccurate extrapolation can yield a reduction in the prediction accuracy.

The extrapolated airfoil is then interpolated based on the relative thickness distribution shown in Figure 3. The interpolation results at various relative thicknesses are shown in Figure 2. The chord distribution of the IEA 15 MW reference wind turbine is applied without any change, but the twist distribution is modified to zero (untwisted blade) to simplify the problem of interest. However, it is understandable that this is not realistic for real wind turbine blades. Therefore, the influence of twist distribution will be given in Section 3.4. The structural properties like mass, bending stiffness, shear stiffness, torsional stiffness and axial stiffness distributions are adopted without any adjustment. This includes the shear center distribution, neutral axis, principal axis orientation and mass axis orientation, but these parameters are not shown in this paper for conciseness.

2.3. Description of the Test Cases

In the present studies, the test case is limited only to cover the blade characteristics to simplify the physical interactions and also to enable isolation of the blade excitation effects. Therefore, the tower and all other structures are set to rigid and only one blade is modeled. The blade is parked at an azimuth angle of zero (blade pointing up) with the blade being pitched to 71.45°, which effectively creates a mean angle of attack of 18.55° with the wind coming from upstream. This mean angle of attack is chosen to match the experimental data for a pitching airfoil case defined in [41]. The rotational speed of the rotor is set to zero to simulate the turbine condition in parked operation. The calculation for the momentum theory is no longer relevant for this case; therefore, the induction calculations and the associated engineering sub-models like tip loss correction and dynamic wake are switched off. However, dynamic stall effects are included since they are highly relevant for the case. In the present studies, three dynamic stall models are tested: Øye, BL and IAG models.

The blade is exposed to sinusoidal variation of the wind velocity in the following form:

$$u\left(t\right)=U_{\infty }cos\alpha \left(t\right)$$

(1)

$$v\left(t\right)=U_{\infty }sin\alpha \left(t\right)$$

(2)

where

$$\alpha \left(t\right)=\Delta \alpha sin\left(2\pi ft\right).$$

(3)

Here, $\Delta \alpha$ represents the excitation amplitude and is set to be constant at 10.45°, which again corresponds to the experimental data for a pitching airfoil case defined in [41]. In the Ohio State University (OSU) experimental campaign [41], various frequencies and two different amplitudes were tested. The amplitude chosen in the present paper is the largest value that will represent different ranges of flow, from attached to fully separated. The values of the incoming wind speed $U_{\infty }$ and excitation frequency f are varied from 10 m/s to 50 m/s and from 0.086 Hz to 1.03 Hz, respectively. Depending on the test cases being investigated, it is often more convenient to represent the frequency in non-dimensional form. In the studies, this is achieved either in form of the reduced frequency:

$$k={\displaystyle \frac{2\pi fc}{2U_{\infty }}}$$

(4)

or as a ratio to the natural frequency. All simulations are performed for 100 s with a timestep size of 0.02 s.

Since the studies are focused on the edgewise resonance, it is a logical choice to use the natural frequency of the blade edgewise mode ($f_0$). For this purpose, a modal analysis was performed in Bladed which yields the modes listed in

Table 1. Results of the modal analysis performed on the blade.

No	Mode Name	Modal Frequency [Hz]
1	1st flapwise mode	0.508
2	1st edgewise mode	0.687
3	2nd flapwise mode	1.515
4	2nd edgewise mode	2.059
5	3rd flapwise mode	3.073
6	3rd edgewise mode	4.035
7	1st torsional mode	4.485

. It can be seen that the 1st edgewise modal frequency is about 0.687 Hz. This means that the frequencies being evaluated (0.086 Hz–1.03 Hz) will be within the range of 0.125$f_0$ to 1.5$f_0$.

3. Results and Discussion

3.1. Comparison of Dynamic Stall Model Against Pitching Airfoil Data

To verify the accuracy of dynamic stall models, Bladed calculations are compared with experimental data of a pitching S809 airfoil. This is performed specifically at the blade station located at 93.68 m. Since the pitching airfoil data obtained from Ref. [41] do not include flexibility (measured on rigid airfoil in a wind tunnel), the flexibility setting for the blade is turned off in Bladed. Equations (1) and (2) are applied through a Wind DLL module in Bladed with an angle of attack variation in Equation (3) being set such that the reduced frequency in Equation (4) at 93.68 m matches the measurement campaign. This effectively simulates a pitching airfoil case at k = 0.079.

The results for all three force components—lift ($C_L$), drag ($C_D$) and pitching moment about a quarter chord ($C_M$)—are presented in Figure 4. Both static and dynamic datasets obtained from the experimental campaign are presented for comparison. It can be seen that the hysteresis effects of $C_L$, $C_D$ and $C_M$ are massively underestimated by the BL and Øye models. The IAG model clearly shows a better match with the experimental data in this deep stall situation. The improved prediction accuracy of the IAG model can be attributed to a more accurate estimation of the leading edge vortex formation [14,15].

The dynamic stall process is initiated by the delay of flow separation, which is captured relatively well by all models. However, when the operating angle of attack is sufficiently high, for some airfoils, an intense leading edge vortex will be formed. This vortex then convects downstream toward the trailing edge. As a consequence of the leading edge vortex formation, pressure decreases significantly on the suction side of the airfoil, causing an increase in lift and drag, while the pitching moment becomes more negative. This intense vortex then breaks down and dynamic stall occurs, creating a significant change in lift, drag and pitching moment magnitudes; see [44,45,46,47,48,49,50] for further details. This stall occurs at about α = 23° for the experimental data in Figure 4, before the magnitudes rise again. The IAG model clearly shows a consistent prediction result against experimental data for the leading edge formation, though the stall at about α = 23° is not captured. In Ref. [14], it was demonstrated that this effect can be modeled by adopting the second-order harmonic terms.

3.2. Dynamic Responses of the Aerodynamic Coefficients

In Section 3.2.1 and Section 3.2.2, the effects of blade flexibility are taken into account. Here, the dynamic characteristics of lift, drag and pitching moments are evaluated against two different conditions. First, the same reduced frequency as applied in Section 3.1 will be adopted in the discussion given in Section 3.2.1, but now with the blade flexibility being activated. Second, Section 3.2.2 discusses the dynamic responses of the aerodynamic forces when the excitation frequency f is set to match the first edgewise modal frequency of the blade; see Table 1.

3.2.1. Exposed to a Constant Reduced Frequency

The simulations at a constant reduced frequency of k = 0.079 (defined at 93.68 m blade station) were performed within $U_{\infty }$ = [10,20,30,40,50] m/s. Here, the magnitude of the wind speed becomes important on top of the excitation frequency because this drives the blade deflection. For the sake of conciseness, only the results for $U_{\infty }$ = 10 m/s, 30 m/s and 50 m/s for lift are presented in Figure 5, Figure 6 and Figure 7, respectively. Results for different dynamic stall models for both rigid and flexible blades are presented. The complete plots for lift, drag and pitching moment are available in Appendix B.1.

For the smallest wind speed of 10 m/s, it can be seen in Figure 5 that blade flexibility has practically almost negligible influences. The difference between rigid and flexible settings remains small. Now, by increasing the wind speed to 30 m/s in Figure 6, it can be seen that the results with a flexible blade start to deviate from the rigid blade setting. This is particularly seen in the upstroke regime (e.g., the higher lift part of the hysteresis). This behavior stems from the fact that the angle of attack increases within this part of the dynamic motion and dynamic solutions are more sensitive to engineering models being applied. The IAG model shows increased vibration characteristics compared to other models. Looking back at Figure 4, it seems that this is how the blade will likely react in reality (note that the same reduced frequency is simulated here).

As for the highest evaluated wind speed of 50 m/s, results start to become much more radical, as seen in Figure 7. It can be seen that the Øye model, BL model and IAG model behave fairly differently. The Øye model shows some increased variations of $C_L$ in terms of the force magnitudes and the angle of attack range. In contrast to the Øye model, the BL model shows somewhat reduced hysteresis curves. A similar behavior is also observed in the IAG model to some degree. However, since the vortex effects in the IAG model are more pronounced than the BL model, it can be seen that force fluctuations are actually increased locally at high angles of attack. The average dynamic polar gradient of the IAG model seems to increase, which might be attributed to the impulsive effects.

3.2.2. Exposed to Edgewise Resonance Frequency

In this section, the blade is exposed to a constant frequency that matches the first edgewise modal frequency of the blade, which excites resonance. The simulations were performed at wind speeds varying from 10 m/s to 50 m/s. This means that the reduced frequency itself decreases with the wind speed. Similar to the discussion in Section 3.2.1, only the results for $U_{\infty }$ = 10 m/s, 30 m/s and 50 m/s for lift are presented for the sake of conciseness. The complete plots for lift, drag and pitching moment are available in Appendix B.2.

In Figure 8, the results for $U_{\infty }$ = 10 m/ are presented for three different dynamic stall models. It shall be noted that this case has the highest reduced frequency of k = 0.598, which in theory shall generate the strongest vibration. In contrast, the wind speed is the smallest (10 m/s), which limits the flexibility effects. It can be observed in Figure 8 that the Øye model yields the weakest vibrations. The IAG model is characterized by the strongest hysteresis effects. This is mainly caused by the impulsive effects being dominant at high reduced frequency. On top of that, the blade flexibility, albeit at the relatively small wind speed, enhances the dynamic stall effects which yield increased vibrations.

Increasing the wind speed to 30 m/s in Figure 9 yields a drastic change of the vibration characteristics. This test case has a reduced frequency of about 0.199 that is still considerably high. On top of that, blade flexibility effects start to become stronger. It can be seen that the vibrations in the Øye model are amplified massively, which occurs mainly only during the vortex breakdown phase; note the discrepancies of the rigid with the flexible case occur mainly starting from α = 18°. The BL model is severely affected by the resonance effects throughout the entire considered excitation range. It can be seen that $C_L$ varies violently from the lowest to the highest angle of attack. The IAG model, on the other hand, shows a somewhat more stabilized vibration behavior and the gradient of the dynamic polar is larger than the other models. At this reduced frequency, it is expected that the strength of the impulsive and circulatory effects is fairly comparable.

For the highest considered wind speed of 50 m/s in Figure 10, the vibrations from the Øye and BL models are now stabilized to a smaller level than the 30 m/s case. This comes from the fact that the actual reduced frequency is smaller than the 30 m/s case compared to the 50 m/s case. The vibration characteristic for the IAG model is fairly comparable with the 30 m/s case. For this flow condition, the circulatory effects become stronger while the impulsive effects become smaller.

One of the longstanding issues that wind turbine designers face is the underestimation of the aerodynamic damping in the aeroelastic modeling of wind turbines. When compared to what the engineers observed in the field, the predicted vibrations are often too exaggerated. It can be seen in the present studies that both BL and Øye models yield massive vibrations for the 30 m/s case, while the IAG model shows a lower estimated value. Since there are no experimental data available for this case, the exact accuracy of the model cannot be assessed. However, the trend different models produce is clearly demonstrated. This outcome may be employed for providing a first quantified measure about the model behavior in the excited edgewise conditions and for selecting a dynamic stall model.

3.3. Angle of Attack Range Variations

In this section, the range of the angle of attack change from the flexible blade simulation to the rigid blade simulation is compared. This is formulated as:

$$\Delta {\alpha }^{Flex}-\Delta {\alpha }^{Rigid}=\left({\alpha }^{Flex,MAX}-{\alpha }^{Flex,MIN}\right)-\left({\alpha }^{Rigid,MAX}-{\alpha }^{Rigid,MIN}\right),$$

(5)

which highlights how much blade vibration grows over the time. Similar to the discussion in Section 3.2, two cases are considered: first, at the same reduced frequency as applied in Section 3.1. For the second case, the excitation frequency is varied from low to high frequency levels within $0.125\le f/f_0\le 1.5$.

3.3.1. Exposed to a Constant Reduced Frequency

Figure 11 displays the vibration growth for the case at a constant reduced frequency of k = 0.079 (calculated at 93.68 m of the blade length). It can be seen that the range of the angle of attack is not altered too significantly for the low-wind-speed case up to 30 m/s. A larger difference between the rigid case and the flexible case starts to occur starting from the wind speed of 40 m/s. It can also be observed that the IAG dynamic stall model has the largest increase in range up to 40 m/s, before it decreases at the wind speed of 50 m/s. On the other hand, the BL model consistently shows the smallest level of range change across all wind speeds. It is also notable that the case without the dynamic stall model shows a consistent increase in the range with increasing wind speed.

3.3.2. Exposed to Varying Frequencies

The aeroelastic responses of the blade section at varying frequencies are discussed in the present section. To enable a clear comparison of the relative relationship between the excitation frequency (f) and the first edgewise natural frequency ($f_0$), the excitation frequency is reported as the normalized frequency ($f/f_0$), and this ranges from 0.125 to 1.5. The results are plotted in two different ways. Figure 12 groups the plots based on the excitation frequency and the data are plotted with respect to the wind speed, while in Figure 13, the results are grouped based on the wind speed and the relationship is plotted against the excitation frequency itself.

Figure 12 clearly illustrates how the reduced frequency affects the vibration characteristics of the blade. Note that the same excitation frequency yields a different reduced frequency (k) for a different wind speed. At the excitation frequency level $f/f_0$ of 0.25 in Figure 12a, one can see that the IAG model generally shows the largest deviation compared to the rigid blade case before it is overtaken by the no dynamic stall case at 50 m/s. Here, it can be seen that the BL model case generally shows the lowest deviation compared to the rigid case. By increasing the excitation frequency to half of the first edgewise natural frequency ($f/f_0$ = 0.25) in Figure 12b, the IAG model now shows the lowest deviation against the rigid case, followed by the BL model, Øye model and no dynamic stall model cases. This observation is mostly pronounced with increasing wind speed up to 40 m/s, before it decreases again when reaching 50 m/s case. The reason is associated with the reduction in the reduced frequency level becoming more significant than the blade flexibility effect itself.

A similar story is found at the non-dimensional excitation frequency of 0.75 in Figure 12c, where the discrepancy against the rigid case increases up to a certain wind speed. Here, it is also observed that while all other cases generally show an increased level of discrepancy against the rigid blade case, the BL model is largely not affected by this particular frequency. Having observed the magnitude of the frequency in more detail, the applied frequency is actually close to the first flapwise natural frequency of the blade. Therefore, the flapwise mode is actually what is being excited.

Now, looking closer at the actual resonance frequency of the first edgewise mode in Figure 12d, the IAG model shows the least deviation compared to the rigid blade case and is followed by no dynamic stall, Øye model and the BL model cases. This frequency excitation is considered the most challenging case because generally, the edgewise mode is prone to negative damping and aeroelastic codes tend to underestimate the aerodynamic damping. The IAG model clearly shows an increased damping level under this condition.

Further increasing the excitation frequency in Figure 12e,f, it can be seen that, generally, the results between the rigid and flexible cases do not deviate too much. This shows that the effects of resonance are not prominent and the blade remains in stable conditions.

In summary, Figure 12 shows that different models have different trends when the wind speed increases, depending on the applied excitation frequency. Except at resonance, generally either the BL model or the IAG model shows the lowest $\Delta \alpha$ increase. This effect is mainly attributed to the flow separation effects in these models. The BL and IAG models use similar flow separation equations, and thus, the results for $\Delta \alpha$ variation are generally similar. However, the actual force hysteresis is actually massively different due to differences in the vortex lift terms (not shown here for the sake of conciseness). Nevertheless, this does not create a massive change in the $\Delta \alpha$ growth. When resonance occurs, the vortex lift effect in the IAG model stabilizes the vibration and actually limits the growth of $\Delta \alpha$. In instability predictions of wind turbine blades, self-excitation of the blades enhances the edgewise modes, and when adopting the Øye and BL models, one can see that strong vibrations might be predicted. In contrast, the vibrations might still occur when using the IAG model, but at lower magnitudes.

Figure 13 is intended specifically to assess the characteristics of the dynamic stall model on the aeroelastic response for varying frequencies. At the smallest wind speed case of 10 m/s in Figure 13a, the excitation frequency largely has no impact on the range of angle of attack. The range only slightly increases for the no dynamic stall case and the IAG model case at the edgewise resonance frequency. By increasing the wind speed to 30 m/s, the blade flexibility plays a stronger role in the aeroelastic characteristics, as shown in Figure 13b. It can be clearly observed that the excitation occurs at $f/f_0$ = 0.5 and $f/f_0$ = 1.0 due to edgewise excitation, while it occurs at $f/f_0$ = 0.75 due to flapwise excitation. For the BL model and Øye model cases, the blade section is excited strongly at the resonance frequency and this leads the angle of attack variation to grow significantly.

Further increasing the wind speed to 50 m/s in Figure 13c, it is notable that most cases deviate from the rigid blade case for all frequency ranges, except for the BL model. However, it is clearly observed that the deviation is smaller than the one observed in Figure 13b. The explanation lies in the fact that the blade flexibility effects are prominent for this wind speed case, but the excitation reduced frequency itself is actually smaller than the case with a wind speed of 30 m/s. This is also supported by the discussion in Figure 12d.

3.4. Twist Angle Distribution Influence

It has been discussed above that dynamic stall model plays a significant role in the aeroelastic effects depending on the excitation frequency. Now, to gain deeper insights into design considerations, it is interesting to see if twist angles affect the characteristics. Therefore, the original twist angle distribution in Figure 3 is now applied. First of all, to evaluate how much twist angle changes the modal properties (keeping all other structural constraints the same), modal analysis of the blade is carried out and listed in

Table 2. Results of the modal analysis performed on the blade for two different twist distributions.

No	Mode Name	Modal Frequency—Zero Twist [Hz]	Modal Frequency—Twisted [Hz]
1	1st flapwise mode	0.508	0.509
2	1st edgewise mode	0.687	0.685
3	2nd flapwise mode	1.515	1.516
4	2nd edgewise mode	2.059	2.054
5	3rd flapwise mode	3.073	3.078
6	3rd edgewise mode	4.035	4.028
7	1st torsional mode	4.485	4.489

. It can be seen that the modal frequencies are altered slightly, but not to a significant degree. The same inflow conditions are then imposed, but the blade set angle is adjusted to maintain the same angle of attack range at r = 93.68 m for the rigid blade case in comparison with the untwisted blade.

The aeroelastic response of the blade is presented in Figure 14 for different dynamic stall modeling approaches. For all cases, it can be seen that the addition of twist reduces the growth of the angle of attack range when the blade flexibility is considered. The effects are mostly prominent for the no dynamic stall case and the case employing the Øye model. This is followed by the IAG model and the weakest effect is observed for the BL model. The twist distribution changes the acting angle of attack and pushes the operating range further away from the stall. The percentage of change is also shown in Figure 14. This shows the relative change of $\Delta \alpha$ when twist is applied. The relative change of the amplitude differs depending on the model. It is noted that the BL model shows the largest absolute value of $\alpha$ growth because the model is massively excited in resonance excitation. This is similar to the results observed in Figure 12d. When applying twist, although the rigid case operates at a lower mean angle of attack, the BL model is still in the excited phase, and thus, the reduction in the $\alpha$ growth is minimal. For the other dynamic stall models (Øye and IAG), the relative reduction is actually at a similar percentage of change (but not in the absolute terms). From the observation, it can be seen that the twist angle distribution will likely become an important aspect to enhance large wind turbine blade stability during parked or idling conditions. When designing the blade, twist distribution is usually aimed at obtaining the required turbine performance by redistributing the angle of attack. However, a more holistic analysis in performing design load cases may be adopted by considering edgewise excitation or known critical instability cases where twist optimization can be performed.

4. Conclusions and Outlook

The present studies have been conducted to evaluate the sensitivity of different dynamic stall models to edgewise vibrations occurring on large flexible wind turbine blades. The investigations were conducted systematically by starting from the rigid response of the blade section with respect to oscillating wind conditions. The blade flexibility effect is then taken into account, causing a dynamic response of the angle of attack. From the discussions, the following aspects can be summarized:

• It is clear that dynamic stall model plays a strong role in the aeroelastic characteristics. The angle of attack range starts to grow away from the rigid blade case when the excitation frequency is in the range of $0.5\le f/f_0\le 1.0$. It is observed that in this range, the flapwise mode is one source of the excitation, though the strongest cause is related to the edgewise mode at $f/f_0$ = 0.5 and $f/f_0$ = 1.0. This can be used as a rough guidance when performing wind turbine design and running load cases, especially in idling or parked conditions.
• The IAG model shows the best agreement with the pitching airfoil data, which represents the rigid conditions. Under the influence of blade flexibility, the IAG model is most sensitive to frequency excitation far from the edgewise resonance. However, when the excitation is close to the edgewise resonance frequency, the IAG model shows a higher aerodynamic damping, which increases the blade stability. This explains the characteristics of this model in generating lower loads for the certain idling cases when self-excitation of the blades enhances the edgewise modes.
• The characteristics of the IAG model for different excitation frequencies are also dependent upon the impulsive contributions. A higher reduced frequency typically increases the impulsive terms.
• Twist distribution is crucial in determining the aeroelastic stability of the blade. Adding a twist distribution (assuming all structural properties are the same) only slightly alters the blade modal frequencies, but affects the operating range of the angle of attack. A more holistic analysis in performing design load cases may be adopted by considering edgewise excitation or known critical instability cases where twist optimization can be performed.

Future studies might be aimed at investigating the dynamic stall sensitivity at higher excitation frequency levels to gain a better picture of the higher harmonic effects. Furthermore, incorporating the torsional effects will also be important to evaluate the coupling between the edgewise and torsional terms and how the aerodynamic damping actually influences the aeroelastic characteristics.