Wind turbine blades are complicated structures that undergo cyclical rotating loads in stochastic wind conditions, which makes wind turbines a complex, non-linear dynamic system. Modeling wind turbines require many degrees of freedom (DOFs) to capture all features of the dynamic behavior of the rotor. The increasing size of the state-of-the-art rotor and these interlinking factors make wind-tunnel studies of next generation super-turbines virtually impossible to extrapolate to the prototype scale, and thus prompt the need for full-scale studies using computer models. These models should also be capable of handling the complexities introduced by adding control systems that are coupled with aeroelastic modes of operation.
The present study uses a novel numerical model capable of handling the aforementioned complexities. We provide a brief description of the main features of the Dynamic Rotor Deformation–Blade Element Momentum (DRD-BEM) model [
19], which is capable of representing the coupled phenomena of wind turbines using two advanced numerical schemes. Structural response is modeled using an innovative approach for heterogeneous composite blades that represents complex modes of blade deformation while reducing the computational efforts substantially [
20]. A dimensional reduction technique for complex beams originally proposed by Prof. Hodges and his collaborators [
21,
22] known as the Generalized Timoshenko Beam Model (GTBM) approach is used in effect with some modifications. Second, the flow behavior is modeled using a novel aerodynamic technique based on the Blade Element Momentum (BEM) approach that allows representation of the instantaneous deformed configuration and the effects of rotor deformation on the computation of aerodynamic loads. This is achieved by the transformation of velocities acting at the rotor level through a series of orthogonal matrices projecting them on to the blade section, and in the same way re-projecting the forces acting at the blade section back to the rotor level where the interference is applied. These numerical schemes work in the context of a multi-physics solver called the Common Ordinary Differential Equation Framework (CODEF) [
19], which also includes modules that model the dynamics of the control system and electromechanic devices on the drive-train.
2.2. Rotor Flow Model: DRD-BEM
The DRD-BEM is a numerical model representing the flow behavior on a wind turbine rotor based on the classical Blade Element Momentum model (BEM) [
25]. Whereas the classical BEM approach misrepresents the aerodynamic forces due to a lack of complete representation of deformed cross-sections of the blade, the DRD-BEM model accounts for these deformations, and the resulting effects on aerodynamic loads itself. This is achieved through a series of orthogonal matrices that transform the velocities, and aerodynamic loads between coordinate systems representing specific aspects of the turbine. These aspects could be a pre-conformed manufacture setting such as a twist, cone angle, or pre-bending, or associated with a control action such as the pitching angle, yaw rotation, or rotor tilt.
We shall start from the velocity vector of the flow passing through an annular actuator disk aligned with the hub coordinate system,
h. The components of this velocity vector are affected by an axial induction factor
a, and a tangential induction factor
. These represent, respectively, the axial velocity deficit and the tangential velocity increase across the actuator disk. Then the velocity vector,
, of the wind going through an annular actuator aligned with the hub coordinate system
h is given by
where
is the undisturbed wind velocity field referred to the hub coordinate system as shown in
Figure 1a,
is the angular velocity of the rotor, and
is the instantaneous radial distance of an annular disk traced by the rotating blade section (identified by
in
Figure 1b) on the rotor plane. This 3-D construction of
reflects how the stream-tubes associated with each blade element, aligned with
, are deflected by the forces exerted on them by the annular actuators. A set of orthogonal matrices transform the wind velocity defined in a coordinate system aligned with the wind itself,
, into the hub coordinate system,
, to account for cases such as rotor tilt, and changes in yaw or wind direction. This unperturbed wind velocity in the hub coordinate system is obtained as:
where
accounts for the misalignment between wind direction and nacelle orientation due to yaw,
accounts for the misalignment due to rotor tilt through a rotation around the horizontal axis of the
nacelle system, and the azimuthal orthogonal matrix
transforms the wind velocity into the
hub coordinate system
h, by a rotation around the main shaft to the blade instantaneous position.
To compute the relative velocity of wind as seen by the instantaneous deformed blade section,
is projected through a few coordinate systems. This deformed configuration (
,
,
) is defined along the deformed reference-line
l as depicted in the right panel of
Figure 1, and the velocity of wind in this reference-frame,
is given as:
In Equation (
3),
denotes the transformation accounting for coning of the rotor, and
performs a rotation around the pitching axis of the blade. These two transformations result in the
blade coordinate system indicated by the subscript
b as per the IEC standards [
26]. For a detailed description of the concept of coning rotors see [
27,
28,
29]. Consequently, the orthogonal matrix
transforms the velocity of wind from blade coordinate system
b to the non-deformed configuration system
L, defined along the original reference line. Furthermore, the orthogonal matrix
transforms from
L to instantaneous deformed configuration system
l, which is obtained from the solution of kinematic aspects of the structural model presented in
Section 2.1. The intrinsic system
L is aligned to the blade section in the chord-normal, chord-wise, and span-wise directions, and represents the longitudinal axis of the beam in its original configuration, as depicted in
Figure 1c. In addition,
denotes the blade section vibrational velocities coming from the structural model, and
is the velocity associated with the motion of the blade section due to combined action of mechanical devices such as yaw, pitch, and azimuthal rotation. Both these velocity vectors are already expressed in the
l system. As an example of all the transformations, the pitching rotation matrix is shown here,
where
is the pitch angle, which defines the angle of misalignment due to blade pitch that is either an operational parameter or the result of a control action.
Once wind velocity is transformed into the instantaneous blade section coordinate
l, aerodynamic lift and drag forces can be computed using the sectional lift and drag coefficients, and the component of wind velocity in the blade section plane. Aerodynamic coefficients are defined for airfoil profiles (representing selected blade sections) based on the relative angle of attack,
. Now, with the knowledge of the magnitude of wind velocity relative to the blade section,
, and its angle of attack
, the sectional lift and drag forces per unit length of span are computed as,
where
and
are the lift and drag coefficients,
is the air density, and
c is the chord length of the airfoil section. These non-dimensional coefficients quantify the aerodynamic behavior presented by the airfoil profiles and primarily depends on the angle of attack,
. When the airfoils are equipped with a Flow-control Device (FCD), an updated set of aerodynamic coefficients are used based on an additional parameter defined by
. More about this will be discussed in
Section 2.4. The total aerodynamic load acting on the sectional blade element aligned with relative wind direction has components corresponding to the lift and drag forces and is given by
where
is the span of the sectional blade element as shown in
Figure 1b.
These forces aligned with direction of wind incidence are then projected onto the chord-normal and chord-wise directions before being projected back to the
h coordinate system. The aerodynamic loads on the blade element expressed in the
h coordinate system is given as
where
is the transformation matrix to project lift and drag forces onto the chord-normal and chord-wise directions, aligned with the coordinates of
l. The above expression (
8) can also be written as
, or in component form as
where
.
A major step in this model is equating the forces obtained from the Blade Element Theory to the change of momentum defined as per the Momentum Theory. Components of
are hence equated to the rate of change of momentum through the corresponding annular actuator. The component normal to the actuator,
is equated to the change in momentum on
, which is associated with the axial interference factor
a (see expression (
1)), and the tangential component
is equated to the corresponding momentum change associated with tangential induction factor
. A set of equations are used to determine these interference factors in an iterative process for each blade section at every time-step of the analysis, adopting an advanced optimization algorithm to improve the stability and the speed of convergence of the iterative process.
The final step in the process of transitioning between the aeroelastic and structural modes is to compute the distributed loads and moments acting on the blade structure per unit span length. These forces expressed in the deformed configuration system
l constitute both aerodynamic forces and gravitational loads and are required as inputs for the structural model. After determination of the induction factors, the process from Equation (
1) through Equation (
3) is repeated to compute the aerodynamic forces on each blade section, however, this time in the
l system. That is,
, whose first two components will give the chord-normal and chord-wise aerodynamic loads. Additionally, the aerodynamic moment on the airfoil section per unit span length acting around the first axis of
l is computed at the blade sections as
where
is the aerodynamic pitch coefficient, which is the third non-dimensional characteristic of airfoils. As with the lift and drag coefficients,
also now depends on both the angle of attack
and the FCD parameter
. The dependence of aerodynamic coefficients on these two parameters can be observed in an example of NACA 64
3-618 presented later in
Figure 2,
Section 2.4.
2.4. FCD Control Module
Fractional devices on turbine blades that can be triggered to alter the airflow dynamics near the rotor are known as flow-control devices (FCD). Use of such devices in controls bring a two-fold advantage of the ability to vary the control parameter for a range of values while making use of minimal power to execute the control action. An FCD essentially alters the aerodynamics of blade sections based on its relative configuration with respect to an airfoil section via the value of a control parameter. For instance, this could be the length of a micro-tab extension, the relative angle for a leading-edge slat, or the angle of actuation of a trailing-edge flap. The FCD module integrated in CODEF is designed to update the aerodynamic properties of each blade section based on the value of the control parameter that defines the relative configuration for the device, which is provided by the control module. In the example selected for this study, we defined the FCD control parameter as the angle of actuation of a trailing-edge slotted flap, .
Fractional-chord trailing-edge flaps that can be fitted as modular attachments onto existing benchmark blades are of primary interest in the current study. The relative positioning and configurations of trailing-edge flaps play a significant role in modifying the flow around the airfoil, and hence the aerodynamic behavior. Among other properties, slotted flaps have the ability to revitalize the boundary layer on the upper airfoil surface to prevent separation near the trailing edge of blades sections if the relative positioning of the flap with respect to the original airfoil section is properly selected.
Figure 4 presents a schematic of the configuration adopted for a modular device externally attached on the turbine blade. This configuration for airfoil–flap assembly was adopted after an extensive study on optimization of the relative positioning of the
-chord Clark Y profile flap near the trailing edge (see Menon et al. [
30]). The relative dimensions and positioning for the flap attachment are shown in
Figure 4a.
Figure 4b shows an example of the implementation of this configuration on a modified NACA 64
3-618, and
Figure 4c presents a blade span view of a modular device externally attached on the turbine blade based on this concept.
Aerodynamic loads acting on the blade is the cumulative effect of forces and moments acting along each section of the blade, which are given by Equations (
5), (
6) and (
10) (see
Section 2.2 for details). At the sectional level, these forces are primarily determined based on the non-dimensional aerodynamic coefficients of lift
, drag
, and pitching moment
. These coefficients characterize each airfoil section, whose values vary with changing angle of attack
. The aerodynamic forces acting at blade sections are modified when airfoil sections are fitted with trailing-edge flaps, and the behavioral alterations depend on the flap chord, flap span, and flow alterations based on the airfoil–flap configurations. Fractional FCDs such as the trailing-edge flap are capable of modifying airflow near the airfoil tail, causing noticeable variations in the aerodynamic characteristics of the airfoil sections and providing a new set of aerodynamic coefficients for each airfoil–flap configuration.
The modified aerodynamic coefficients were computed for two key airfoils typically used in the outer regions of the blade span: the NACA 64
3-618 and the DU 93-W-210, when attached with a
-chord Clark Y profile trailing-edge flap. These airfoil sections are among the more aerodynamically efficient sections and are widely used in contemporary wind turbine blade designs, such as the benchmark wind turbine designed by National Renewable Energy Laboratory (NREL), known as the NREL-5MW Reference Wind Turbine (RWT). On such a blade, these two airfoil sections cumulatively make up about
of the span, as indicated in
Figure 5. The inner regions of the blade (closer to the root) have airfoils that are thicker to ensure structural stability, whereas the outer regions (closer to the tip) use thinner airfoils, which contribute more to the aerodynamic efficiency of the blade when compared to the inner regions. A major share of this aerodynamic contribution to blade operation originate in the shaded regions in
Figure 5, which are essentially the span region equipped to be attached with a trailing-edge flap. The relative positioning of the trailing-edge flap adopted in this study is depicted in
Figure 4. As mentioned earlier, the configuration of the airfoil–flap assembly plays a key role in determining the quantitative modification in aerodynamic behavior. These airfoil sections were studied for a range of configurations of the airfoil–flap assembly, and is defined using the relative angle
between the airfoil and flap chords. The repository for aerodynamic characteristics of these airfoil section are available for a range of
to
, evaluated at regular intervals of
.
The control system module in CODEF currently has the ability to integrate the dynamics of control techniques such as yaw, pitch, and coning. Using trailing-edge flaps as a prototype, the module is extended with the capability to incorporate the dynamics of flow-control devices (FCD), simulating the interaction of such control actions with the dynamic aeroelastic response of the rotor. This means that the effects on the rotor dynamics due to a control decision from the FCD module and vice-versa will be evaluated at every instant of the simulation. In the example analyzed here, the aerodynamic characteristics of the modified airfoil sections fitted with fractional trailing-edge flaps are made available to the control module, which has the capacity to interpolate the instantaneous values of the aerodynamic coefficients from the repository of aerodynamic data. At each instant of the time-step analysis (see
Section 2.2), the adaptive algorithm evaluates the level of actuation of the trailing-edge flap
for each section of the blade, based on the instantaneous input from the control module. Then, each blade section adopts an updated value for its aerodynamic coefficients. The flowchart shown in
Figure 6 gives an overview of functional algorithm that is used by the control system module of CODEF.
The aerodynamic coefficients and the resulting loads and moments acting at each airfoil section will now depend on two instantaneous parameters: angle of attack at the airfoil section, and angle of flap actuation defining the airfoil–flap configuration. The adaptive ODE framework ensures that structural deformations and their effects on the aerodynamic loads that arise as a result of such aerodynamic alterations, are also considered through the natural integration that CODEF makes of the multi-physics dynamics of the machine. One other aspect to note here is that this technique could be extended to other types of FCDs, besides flaps, just by selecting a different repository of aerodynamic data which provides the modified aerodynamic coefficients in function of the angle of attack and the level of actuation of the FCD in question.
Figure 2 shows the non-dimensional sectional coefficients of lift, drag, and pitching moment characterized for NACA 64
3-618. These numerical computations and the resultant aerodynamic coefficients of two key airfoil sections with attached trailing-edge flaps present valuable data for the design of future innovative turbine blades with active FCDs. These characteristics were obtained from an extensive study on two-dimensional flow characterization for these airfoil–flap assemblies. A steady-state pressure-based computational fluid dynamic solver was used to this effect, and the range of flap configuration covers a substantial set of scenarios relevant to wind turbine operating conditions. For more details of this study, the reader is referred to Menon et al. [
30] and the references therein.